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ABSTRACT 

The Ormsby model of h\man dynamic orientation, a discrete 
time computer program, has been used to provide a vestibular 
explanation for observed differences between two washout schemes. 
These washout schemes, a linear washout and a nonlinear wash- 
out, were subjectively evaluated by Parrish and Martin, They 
found that the linear washout presented false rate cues, caus- 
ing pilots to rate the simulation fidelity of the linear scheme 
much lower than the nonlinear scheme. By inputting the motion 
histories from the Parrish and Martin study into the Ormsby mo- 
del, it was shown that the linear filter causes discontinuities 
in the pilot's perceived angular velocity, resulting in the sen- 
sation of an anomalous rate cue. This phenomenon does not oc- 
cur with the use of the nonlinear filter. 

In addition, the suitability of the Ormsby model as a sim- 
ulator design tool was investigated. It was found to be a use- 
ful tool in predicting behavior of simulator motion bases, even 
when the mechanical motion base is replaced by a computer sim- 
ulation, Further investigation of the model could provide -sim- 
ulation designers with a tool to predict the behavior of motion 
bases still in the drawing board stage. 
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CHAPTER I 
INTRODUCTION 


- For many applications it is often desirable to simulate 
a particular vehicle motion without using the actual vehicle 

• The Federal Highway Department sponsors many 
drunk driver studies. In order to insure the 
safety of the driver# the vehicle and the ex- 
perimenters, these experiments are often car- 
ried out in a moving base simulation of an 
automobile. 

« The U.S. Navy has commissioned studies of the 
habitability of large high-speed surface-eff— 
ect-ships. It is necessary to understand to 
what extent crews will be able to fiinction on 
these ships even before a prototype is built. 

This research is carried out on a motion gen- 
erator, which simulates the expected range of 
motion of these ships [7] . 

• The U.S. Air Force makes extensive use of both 
stationary and moving base aircraft simulators 
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in pilot training programs. Simulators pre- 
sent no risk to the pilot, and avoid the costs 
of fuel and repair or possible loss of an air- 
craft . 

The above examples illustrate three of the many possible 
uses of simulators - to carry out driver-vehicle studies with- 
out using an actual vehicle, to predict crew habitability on 
board a ship not yet built, and to train aircraft pilots with- 
out risking the pilot or the plane. As vehicles become in- 
creasingly complicated, and costs continue to rise, motion 
simulation takes on a new importance. 

There are many types of cues a person uses to sense motion. 
The basic inputs are specific force and angular acceleration, 
which can influence the vestibular system in the inner ear, the 
tactile sensors at points of contact with the vehicle, and the 
proprioceptive sensors as muscles are -stretched and compressed. 
In a simulator, it isi not always possible to reproduce a par- 
ticular motion history exactly. Often, some cues can be pre- 
sented only at the expense of neglecting other cues. The basic 
goal in motion simulation is to arrive at a compromise in pre- 
senting the cues, in order to best represent the desired motion. 
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1«1 The Physiology of Motion Simulation 

Sinmlation technology now makes heavy use of digital 
computers to present as much of the motion cue as possible. 

High speed processing allows the use of very complex linear 
filters, and recently, of nonlinear adaptive filters. Micro- 
processor technology has also made much of the slower elec- 
trical circuitry obsolete. 

But the goal of simulation has not really changed - try 
to present as many of the specific force and angular acceler- 
ation cues as possible, without exceeding the constraints of 
the simulator [18] . This has always been the most straight- 
forward approach, since it is the specific force and angular 
acceleration cues which are most readily available. 

Once a good understanding of the physiological aspects 
of motion simulation is attained, a physiological model of the 
human operator will be a valuable tool in simulator design. 

The comparison of actual motion and simulated motion using such 
a model would be useful in determining'' the realism of the sim- 
ulation in a quantitative way. This model would also be help- 
ful in comparing two different simulation schemes, providing 
a quantitative measure of their differences. 
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1.2 The Use of Washout and Visual Cues in Simulation 

Constraints in position, velocity and acceleration of a 
simulator limit the capability of producing a desired motion 
exactly. The problem is to present the sensations of a wide 
range of motion, and to do this in a very limited space. This 
problem is solved with the use of washout filters in each axis 
of motion, in order to attenuate the desired motion until it 
falls within the constraints of the simulator. 

An important aspect of motion simulation has not yet been 
mentioned - the visual cues available to detect motion. Peri- 
pheral visual cues seem to be most important in presenting the 
sensation of motion. The peripheral field may be stimulated 
by a moving pattern of stripes or dots, or by an actual "out - 
the - window" cockpit view [2,5], 

Taken together, washout filters and visual stimulation 
perform the function of simulation in which motions seem to go 
beyond the constraints of the simulator. The motion is dupli- 
cated to the point of constraint in a given axis. Then the wash- 
out filter takes over and attenuates the motion to meet the 
constraint. Meanwhile, the visual field is stimulated so as to 
give the impression of continued motion, motion beyond the cap- 
abilities of the simulator. In this way, a wide range of mo- 
tions can be simulated using a very restricted motion base. 
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1«3 Thesis Objectives and Organization 

It is obvious from the previous discussion that the wash- 
out filters in a simulator are critical to the fidelity of the 
simulation. The research leading to this thesis compares two 
different types of washout filters currently in use, in order 
to quantify the differences between them,. The means of com- 
parison is a physiological model of human dynamic orientation, 
based largely on the known physiology of the vestibular system. 
This work attempts to answer a specific question and a general 
question: 

• Can the observed differences in simulation 
fidelity between the two filters be explain- 
ed using a physiological model of human dy- 
namic orientation? 

•What are the implications~for this model as 
a drawing board tool in simulator design? 

Chapter II presents the two washout filters in detail, and 
discusses the previous work which led to the research present- 
ed in this thesis. 

Chapter III describes the human vestibular system and the 
model of human dynamic orientation developed by Ormsby. 

Chapter IV describes the data in this work, as input to 
the model, and then presents the perceived angular velocities 
as output from the model. 
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Finally, Chapter V presents the conclusions which can be 
drawn from the results presented in Chapter IV, in light of 
the questions posed in the above- section. Also included are 
suggestions for further research in this area. 
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CHAPTER II 
THE WASHOUT FILTERS 

The two washout filters of interest in this comparitive ‘ 
study are the following: 

• A linear filter, essentially a Schmidt -and' Conrad 
coordinated washout [16,17], 

• A nonlinear filter, coordinated adaptive washout. 

Basically, the two filters are versions of Schmidt and Conrad's 
coordinated washout. This scheme uses washout filters in the 
three translational axes, and only indirectly washes out the 
angular motion. The primary difference between the linear and 
nonlinear schemes is in the type of translational washout fil- 
ters employed. The linear scheme uses second-order classical 
washout filters in the three axes, while the nonlinear scheme 
uses coordinated adaptive filters for longitudinal and lateral 
washout and digital controllers for vertical washout. These 
schemes differ in their presentation of the rate cues, for a 
pulse input. The linear scheme presents an anomalous rate cue 
when the pulse returns to zero. This behavior is not observed 
with the nonlinear scheme. 
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The next two sections discuss the filters in greater de- 
tail, The final sections present the differences between the 
filters and the results of a previous subjective analysis of 
the washout schemes. 

2,1 The Linear Washout 

-The purpose of washout circuitry is to present transla- 
tional accelerations and rotational rates of the simulated air- 
craft. It is necessary to obtain coordination between trans- 
lational and rotational cues in order to accomplish certain 
motion simulations: 

• A sustained horizontal translational cue can 
be represented by tilting the pilot. The 
gravity vector is then used to present the 
cue. But in order to make this process be- 
lievable, the rotation necessary to obtain- 
the tilt angle must be below the pilot’s ab- 
ility to perceive rotation. The solution is 
to start the cue with actual translational 
motion of the simulator until the necessary 
tilt angle is obtained. In this manner, the 
pilot will sense only translational motion, 
long after such motion has actually ceased. 

• In a similar sense, it can be seen that a de- 
sired roll or pitch cue cannot be represented 


22 


by means of rotation alone. This would result 
in a false translational cue, because the gra- 
vity vector is misaligned. In order to present 
a rotational cue, translational motion must be 
used at the start, to offset the false trans- 
lational motion cue induced by the rotation. 

The two cases above clearly illustrate the need for coor- 
dination in translational and rotational motion, Schmidt and 
Conrad's coordinated washout scheme fulfills this' need. Fig- 
ure 2.1 presents a block diagram illustrating the basic con- 
cepts . 

The desired motions of the simulated aircraft are trans- 
formed from the center of gravity of the aircraft to the cen- 
troid of the motion base. This transformation provides the de- 
sired motion at the pilot's seat. The motions of the base are 
based on the desired motions of the centroid. 

Vertical specific force is transformed to vertical accel- 
eration 2 ^ by use of a second-order classical washout filter. 

The longitudinal and lateral accelerations are also obtained 
from the longitudinal and lateral specific forces. First, these 
specific forces are separated into steady-state and' transient 
parts. The steady-state part of the cue is obtained from a 
tilt angle to align the gravity vector. The transient part of 
the cue is transformed into the longitudinal acceleration, 
and the lateral acceleration, y^, by a second-order classical 
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washout filter. 

Braking acceleration is then used to keep the motion with- 
in the prescribed position, velocity and acceleration limits 
of the motion base. 

The rotational degrees of freedom are only indirectly 
washed out through elimination of false g cues. Rotational 
rate cues are represented by angular and translational motion, 
just as longitudinal- or lateral cues. But in this case, the 
translational motion is used only to eliminate the false g cue 
induced by rotational movement, and thereby makes no direct 
contribution to the rotational cue. 

After the six position commands (x^,y^,z^,ij; , 0 ,<^) are ob- 
tained from the washout circuitry, lead compensation is pro- 
vided to compensate for servo lag of the base. The actuator 
extension transformation is then used to obtain the correct 
actuator lengths used to drive the motion base. 

The actual filter evaluated in this work is a Schmidt and 

£ 

Conrad coordinated washout, adapted by Langley Research Center 
[14] , The major difference is that the Langley washout is car- 
ried out in the inertial reference frame, rather than the body 
axis system, A block diagram of this filter is shown in Fig- 


ure 2.2 
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Figure 2.2 Block diagreim for Langley linear washout scheme [13] 
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2.2 - The Nonlinear Washout 

The nonlinear filter of interest here is again essentially 
a Schmidt and Conrad coordinated washout. The difference be-' 
tween the nonlinear Langley filter and the Schmidt and Conrad 
filter are that the Langley filter uses the inertial reference 
frame rather than the body axis system, and nonlinear filters 
are used for the washout rather than the linear filters used 
by Schmidt and Conrad. Hence, the designation “nonlinear wash- 
out” is used. 

Figure 2.3 presents a block diagram for this nonlinear 
scheme. It is seen that two different types of nonlinear fil- 
ters are used - coordinated adaptive filters for longitudinal 
and lateral cues, and digital controllers for vertical cues. 
These two types of filters will be discussed in turn. 

Coordinated adaptive filters [11] are based on the prin- 
ciple of continuous steepest descent. They are used in this 
washout scheme to coordinate surge and pitch in presenting the 
longitudinal cues, and sway and roll in presenting the later- 
al cues. Derivation of these filters can be found in the liter- 
ature [11,12] . Basically, they perform the same functions as 
the second-order classical filters used by Schmidt and Conrad 
by providing translational specific force cues and rotational 
rate cues. 

Digital controllers, the second type of nonlinear filters, 
are used to provide the uncoordinated heave and yaw cues. A 
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Figure 2,3 Block diagram for Langley nonlinear washout scheme [13] 
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first-order digital controller provides the yaw rate cue, while 
a second-order controller provides the vertical specific force 
cue. These filters are designed to present as much of the on- 
set cue as possible before switching to the washout logic. 

Figure 2.4 illustrates the design concept for a first- 
order digital controller. From 0 to the controller presents 
a scaled version of the commanded input. At a linear decay 
is applied to reduce; the command to the motion base constraint 
value, B. Washout then occurs at the constrained value, unless 
another input is commanded, as at T 2 * 

The second-order digital controller used for the vertical 
specific force is similar, although mathematically more 
complex.' 


2.3 A Comparison of Washout Schemes 

Essentially, the two washout schemes of interest are 
Schmidt and Conrad washouts. The so-called linear washout 
contains second-order classical washout filters which trans- 
form the specific forces in each axis to translational accel- 
erations in each axis. The Langley washout perfoimis these 
transformations in the inertial frame rather than the body 
axis frame used by Schmidt and Conrad. 

The nonlinear washout scheme uses two types of nonlinear 
filters to provide the translational acceleration cues, A 
coordinated adaptive filter is used to coordinate surge and 




30 


pitch for longitudinal cues, and sway and roll for lateral 
cues. A digital controller is used for the uncoordinated 
heave and yaw motions. Again, the Langley nonlinear scheme 
washes out in the inertial frame. 

In Figure 2.5, amplitude and phase versus frequency is 
shown for the three types of washout filters - linear, adaptive 
and digital controller. Both the first-order and second-order 
cases are shown. The motion base characteristics are the same 
in all cases. Since the amplitude and phase response of the 
nonlinear adaptive filter changes with the magnitude of the in- 
put, the worst case for the nonlinear filter is presented here. 
As is shown, the digital controller has the best response char- 
acteristics, and the adaptive filter is better than the linear 
filter. This holds true for both the first- and second-order 
cases. 

In terras of motion cues, there is a fundamental difference 
between the linear filter and nonlinear filter for the first- 
order case. Figure 2,6 shows the response of the two filters 
to a pulse input. The difference between the filters is the 
anomalous rate cue presented by the linear filter as the pulse 
input returns to zero. This false cue is most noticeable for 
pulse-type inputs, and disappears as the input becomes sinu- 
soidal, Since the differences between the linear and nonlinear 

t • 

filters vary with input, performance of a given filter is depen 
dent on pilot input and simulator responsiveness in each axis. 
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Figure 2.5 Amplitude and phase plots for three types 
of filters [1] 
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Figure 2 . 6 First-order linear and nonlinear adaptive 
response to a pulse input [13], 
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2>4 Empirical Comparison of Washout Filters 

Parrish and Martin, the major investigators of these two 
washout schemes at I*angley, devised a subjective test to deter- 
mine the differences between the two filters in actual simula- 
tion [13] . Seven pilots flew a six-degree-of-freedom simulator 
equipped with both linear and nonlinear washout schemes. .The 
pilots were asked to rate the motion cues presented by each 
scheme for throttle, column, wheel and pedal inputs about a 
straight-and-level condition during a landing approach* 

The results of this evaluation process are presented in 
Table 2,1. Each pilot determined his own criteria for evalua- 
tion. In addition to rating the cues for each input, the pilots 
were asked to rate the overall airplane feel - that is, how 
successful the overall motion was in representing the actual 
airplane. In the table, the open symbols represent the rating 
of the linear method, while the solid symbols represent the 
rating of the nonlinear method. The washout methods were ap- 
plied to a 737 CTOL aircraft simulation, and four of the pilots 
(represented by the triangular symbols) had previous 737 cock- 
pit experience. 

The pilot ratings for the throttle input are the same for 
each method, as shown in Table 2,1, Even given the methods 
back to back for comparison, the pilots could not detect that 
a change had been made. Figure 2,7 shows the time histories 
for such a change in throttle setting. Longitudinal accelera- 
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tion and pitch rate are the inputs ' to the washouts from the 
simulated aircraft for such a maneuver. The figure shows very 
little difference between the washout schemes, as the pilot 
ratings indicated. The fundamental difference between the two 
pitch rate filters is obscured in order to correctly represent 
the decrease in longitudinal acceleration at six seconds. 

An elevator doublet was input to rate the motion cues for 
a column input. Again, the pilots found little difference be- 
tween the linear and nonlinear washout schemes, as shown in 
Table 2,1, Pour pilots rated the filters the same, while the 
other three rated the nonlinear filter slightly higher. The 
time histories for the elevator inputs are shown in Figure 2,8, 
As in the throttle input case, the fundamental difference be- 
tween the pitch rate filters is not apparent, due to the coor-, 
dination between pitch rate and longitudinal acceleration. In 
addition, the pitch response of the 737 is not at all pulse- 
like, which lessens the difference in performance of the filters 

Wheel inputs were evaluated using ailerons to bank the 
simulator 20*^ for a 30° heading change with a return to straight 
and-level flight. The pilots preferred to separate the wheel 
inputs into roll cues and yaw cues to evaluate these cues in- 
dividually. Figure 2,9 shows the time histories for roll cues 
in the maneuver described. The anomalous rate cue is present 
for the linear washout. This is reflected in the pilots* rat- 
ing, as seen in Table 2,1, All seven pilots felt the nonlinear 
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filter to be at least one and one-half categories higher 'than 
the linear filter. 

Figure 2,10 shows the time histories for yaw cues during 
the same aileron maneuver. Again, the anomalous rate cue is 
present for the linear filter scheme. The pilots were parti- 
cularly aware of a negative rate cue when the simulated air- 
craft rate returned to zero during maneuvers of this type. 

The ratings in Table 2,1 are at least one category' higher for 
the nonlinear scheme, reflecting the unnateral feel of the 
linear rate cue. 

Each pilot flew a set of rudder maneuvers for both wash- 
outs to evaluate roll and yaw cues. There were no changes in 
the ratings from those obtained using the wheel. This is re- 
flected in the time histories for roll and yaw, shown in Fig- 
ures 2,11 and 2.12, respectively. 

Finally, each pilot was asked to rate the two washout 
schemes in terms of overall airplane feel. Table 2.1 shows 
the large contribution made by roll represent ationin the over- 
all airplane simulation. All pilots rated the nonlinear wash- 
out at least one and one-half categories higher than the lin- 
ear washout. They specifically objected to the • anomalous rate 
cue presented by the linear filter in both roll and yaw. 

From this study, Parrish and Martin concluded that the non- 
linear washout scheme better represents actual airplane motions 
than does the linear washout method, at least in an empirical 
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sense.. It appears that the nonlinear scheme does not present 
more of the motion cue;. it merely eliminates the false cue pre- 
sent in the use of the linear washout. 

The work presented in this paper attempts to quantify the 
results obtained in the subjective analysis made by Parrish 
and Martin. In order to accomplish this, the motion histories 
from the Parrish and Martin' study are input to a model of human 
dynamic orientation. The output from the model will provide 
a vestibular explanation for the sensation differences between 
the two filters. Results of this work are presented in Chap- 
ter IV. 
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CHAPTER III 

THE PHYSIOLOGICAL MODEL 

A model which . predicts human perceptual response to mo- 
tion stimuli has been developed at M.I.T.'s Man-Vehicle Labor- 
atory by Osrmsby [10] . The model, which exists as a FORTRAN com- 
puter program, is based on the known physiology of the vesti- 
bular system. While little is known about the processing of 
the specific forces and angular accelerations received from the 
vestibular organs, the simplifying assximptions made about this 
process produce a model which agrees with available neurologi- 

X 

cal and physiological data. 

This chapter first presents an overview of the vestibular 
system, and then goes on to discuss the mathematical modelling 
of the system which leads to the current FORTRAN model. More 
detailed descriptions of the vestibular system may be found in 
the literature [9,15,19,20]. The complete derivation of the 
model of human dynamic orientation is found in Ormsby, And a 
description of the actual FORTRAN programs and their use is 
available in the appendix to this thesis. 
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3 ♦! The Hijman Vestibular System 

The vestibular system, or labyrinth, comprises the non- 
auditory portion of the inner ear. It is composed of three 
semicircular canals, one utricle and one saccule in each ear. 

The semicircular canals are the rotational motion sensors. 

They consist of three approximately orthogonal circular tor- 
oidal canals. The canals are filled with a water-like fluid 
called endolyraph. When the head undergoes angular accelera- 
tion, the endolymph tends to lag behind the motion of the canal 
walls. The motion of the endolymph relative to the canal walls 
displaces the cupula, a gelatinous mass which completely ob- 
structs one section of the canal called the ampulla. Sensory 
hair cells embedded at the base of the cupula detect its dis- 
placement, As a result, the deformation of the cupula is trans- 
formed into an afferent firing rate which provides a signal of 
rotational motion to the central nervous system. {see Figure 
3.1) . 

In a particular canal, all of the hair cells have the same 
polarization. When the flow of endolymph displaces the cupula 
in a single direction, the hair cells are either all excited 
or all inhibited. As shown in Figure 3.2, the canals on either 
side are essentially coplanar with the other side. Thus, they 
are pairwise sensitive to angular accelerations about the same 
axis. Since a pair of canals which are sensitive about the 
same axis have opposite polarities, it is assumed that the high- 
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Figure 3.2 Orientation of semicircular canals [10] 
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er processing centers respond to the difference in afferent 
firing rates. 

Two otolith organs, consisting of a utricle and a saccule, 
are located in each ear. The otolith is sensitive to changes 
in specific force. Figure 3.3 depicts the basic structure of 
the otolith organs. The otolith consists of a gelatinous 
layer containing calcium carbonate crystals, known as otoconia. 
This layer is supported by a bed of sensory hair cells. An 
acceleration of the head shifts the otoconia relative to the 
surrounding endolymph, due to the higher specific gravity of 
the otoconia. This shifting causes the sensory hair cells to 
bend, sending a change in afferent firing rate through the af- 
ferent nerve fibers to the central nervous system. 

As shown in Figure 3.4, the utricles are oriented such 
that their sensitivity is in a plane parallel to the plane of 
the horizontal semicircular canals. The sensitivity of the 
saccules is in a plane perpendicular to the horizontal canals. 
The hair cells in the utricle are sensitive in all directions 
parallel to its plane of orientation, while the hair cells in 
the saccule make it predominantly sensitive to accelerations ■ 
perpendicular to the utricular plane. 
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3,2 The Ormsby Model 

The mathematical model of the semicircular canals consists 
of several parts* The first* part is the mechanical model of 
the cupula deflection caused by motion of the endo lymph. The 
second part includes the interaction between the mechanical 
movement and the afferent firing rate. The third part concerns 
measurement noise, which is that portion of the afferent sig- 
nal found to be independent of the- mechanical stimulus input. 

Figure 3,5 depicts the afferent model of the semicircular 
canals as arrived at by Ormsby, Observation of cupula motion 
led to the torsion pendulum model [93V suggesting that the over- 
damped system reacts to angular velocity rather than angular 
acceleration. The results of the model are expressed as a 
transfer function of the following form:- 

FR (s) = (57.3) {300s^) (,01s+l) 

° w(s) 

(ISs+l) (..005S+1) (30s+l) 

+ + n(t) (3.1) 

s 

The model of the otolith system is composed of two parts 
- the mechanical model of the otolith sensor, and the affer- 
ent response to otolith displacement. Figure 3,6 presents the 
afferent model of the otolith system used by Ormsby, The me- 
chanical model of the otolith is that of a fluid- immersed mass 
retained by a spring. The resulting transfer function relat- 
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to(s) - stimulus, angular velocity of head 
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Mechanical coupling of endolyraph 
to effective hair cell displacement 


F = Constant which relates hair cell 
displacement to change in firing 
rate 

SFR ~ Spontaneous firing rate of typical 
afferent cell (90 ips) 


Figure 3,5 Afferent model of semicircular canals (after Ormsby 1101) 







Figure 3,6 Afferent model of otoliths (after Ormsby 110]) 
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ing afferent firing rate to specific force is; 

FR (s) = (18000) (s+.l) 

® SP(s) + SFR + n(t) (3.2) 

(S+.2) (s+200) ® 

The input to the model consists of a stimulus composed of 
specific forces and angular accelerations in each axis of the 
head coordinate system. Each of these afferent inputs is then 
transformed, into sensor coordinates. From this sensor stimula- 
tion, the afferent firing rates are derived, using the trans- 
fer functions presented above. 

At this point, the process becomes purely guesswork. Even 
assuming that these afferent firing rates are available to some 
central processing system in the brain, the form which this 
processing takes is simply a guess. Ormsby guessed that the 
central processor performs a type of least mean squares error 
optimization to make an estimate of the specific force and 
angular velocity inputs based on the afferent firing rates out- 
put from the vestibular system sensors. 

In this case, such a least mean squares estimator is a 
Kalman filter [4,8]. The input is unknown except for an ex- 
pected range of magnitude and a frequency bandwidth, and an ex- 
pected measurement noise. Also, the input and the noise stat- 
istics are time invariant, which makes the filter a steady- 
state Kalman (or Wiener) filter. This steady-state Kalman fil- 
ter is used by the model to produce estimates of specific force 
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and angular velocity from the afferent firing rates. These 
estimates are tuned, using the Kalman filter gains, to yield 
estimates which fit the available neurological and physiologi- 
cal data for known inputs. 

The filters used for canal processing are tuned such that 
the estimates produced for the angular velocities are essential 
ly unchanged from the afferent inputs. This observation is in 
agreement with available data, suggesting that very little 
central processing is performed. The otolith filters must be 
tuned so that a more dramatic effect by the filters on the aff- 
erent input is observed. This suggests that more central 
processing is required, or that the model of the afferent re- 
sponse is missing a term which has subsequently been attribut- 
ed to the central processing mechanism in the tuning procedure. 
Basically, the filter acts as a low pass filter with a time 
constant of 0.7 seconds. The utricle and saccule differ only 
in the Kalman filter gains, where the saccule gains are twice 
the utricle gains. 

Once the specific force and angular velocity estimates 
have been obtained from the Kalman filters, the saccule non- 
linearity must be accounted for. This is done by means of a 
nonlinear input-output function, and allows the model to in- 
clude observed attitude perception inaccuracies known as Au- 
bert or Mueller effects [6] , The resulting specific force and 
angular velocity estimates are transformed back to head coor- 
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dinates . 

These estimates must now be combined to yield new estimates 
of perceived position, velocity and acceleration. In the model 
this is accomplished by a separate scheme, known as DOWN . 
is a vector of length 1 g in the direction of perceived ver- 
tical; as such, it is the model *s prediction of the perceiv- 
ed vertical. The basic assumptions used in combining the spe- 
cific force and angular velocity estimates to arrive at DOWN 
are the following: 

• The system will rely on the low frequency por- 
tion of the specific force estimates provided 
by the otoliths . ' 

• The system will use that part of the canal in- 
formation which is in agreement with the high 
frequency content of the rotational informa- 
tion provided by the otoliths. 

This logic is presented in Figure 3,7, Block A produces 
the estimate of rotational rate from the input specific forces 
assximing ^ is fixed in space. The low frequency component of 
this estimate is filtered out in Block B, Block C isolates 
the component of the low frequency angular velocity estimate 
which is perpendicular to ^ and DOWN . This is the mechanism 
discussed in Chapter II, which allows cancellation of canal 
signals arising when prolonged rotations are stopped sudden- 
ly, The effect of the three blocks is to produce a rotational 
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Figure 3,7 DOWN estimator 
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vector which represents the low frequency rotational informa- 
tion available from the otoliths (R . ) . 

oto 

Block D confirms whether or not the high frequency portion 
of the canal rotational information is consistent with the high 
frequency portion of the otolith rotational information. The 
inconsistent part of the canal information is sent through a 
high pass filter (Block E) and is then combined with the con- . 
sistent portion of the canal inf oirraation . The component of 
the resulting rotation vector parallel to DOWN is then elimi- 
nated, leaving a rotational vector due to canal information 
(R___) . The total estimate of the rotation rate of the outside 

s s w 

world with respect to the last estimate of DOWN , is com- 
puted by subtracting R^^^ from The net result of Blocks 

H and I is to produce an estimate of DOWN which is the same as 
the estimated specific force vector. This is accomplished by 
a slow reduction in the discrepancy between ^ and DOT-TN , elim- 
inating any accumulated errors resulting from the integration 
of rate information. 

Figure 3,8 illustrates the model for predicting perceived 
rotational rate. The angular velocity vector parallel to DOWN 
becomes the perceived parallel angular velocity. The perpen- 
dicular angular velocity is computed in three steps: 

^ / 

1. Calculate the difference between the com- 
ponent of angular velocity perpendicular 
to DOWN, and the angular velocity cons is- 
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tent with the rate of change of the direc~ 

. tion of DOWN (Block K) . 

2, High pass filter this difference, 

. 3, Combine the filtered result with the DOWN- 
consistent angular velocity. 

This process assures that the canals provide the high frequency 

} 

component of the rotational rate, while the low frequency com- 
ponent is the rotational rate consistent with DOWN . The total 
sense of rotation is thus the stim of the parallel and perpen- 
dicular components . 

This completes the description of the form of the Ormsby 
model used in this work. A complete description of the model 
may be found in Ormsby 's thesis. Figure 3.9 presents an over- 
view of the entire model. At this point, a few important ob- 
servations should be made: 

• The Ormsby model was tuned using inputs with 
known outputs for a certain set of discrete 
time intervals - namely, an afferent update 
interval of 0,1 seconds and a Kalman filter 
estimate update interval of 1,0 seconds. In 
this thesis, due to the characteristics of the 
input data, the afferent update interval is 
0,03125 seconds, and the Kalman filter esti- 
mate update interval is 0,25 seconds i In order 
to use these two intervals, the model had to 
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be retuned by changing the Kalman gains. This 
process, which is- necessary each time the update 
intervals are changed, is described in more de- 
tail in the appendix. 

• One important assumption made by this model is 
that the inputs are unknown prior to their pro- 
cessing. It was noted in the introduction to 
this thesis. that specific force and angular ac- 
celeration act on the body as a whole, provid- 
ing visual, tactile and proprioceptive, as well 
as vestibular, cues. This model takes account 
of the vestibular cues only^ although the tun- 
ing process may force it to consider certain 
aspects of the other sensory cues. Thus, when 
this model is applied to cases where the sub- 
ject might have prior knowledge, or at least 
an expectation of the motion, the results must 
be interpreted in light of the limitations im- 
posed by the model. 

The Ormsby model of human dynamic orientation was used in 
this work as a FORTRAN program implemented on a PDF 11/34. 

'The main program, as well as all associated subroutines, is 
documented in the appendix. 
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CHAPTER IV 
DATA AND RESULTS 


As a logical consequence of the two previous chapters, it 
is desirable now to evaluate the two washout schemes using the 
Ormsby model of human dynamic orientation. Such an evaluation 
could serve the purpose of quantifying the differences between 
the two filters which Parrish and Martin found in their subjec- 
tive study. In addition, this evaluation could shed some light 
on the question of the model’s usefulness in simulator design. 

This chapter presents the data used for this study, and 
the results of the processing of the data by the Ormsby model, 

4,1 Data Description 

The data used in this work consists of four runs made 
with a linear or a nonlinear washout -on the Langley simulator. 
These runs coincide with Figures 2,9, 2,10, 2,11, and 2,12, 
Table 4,1 lists the definitions of the variables measured dur- 
ing these simulation runs. Note that not only are the simula- 
tor motions recorded, but also the commanded motions of the 
aircraft. This allows evaluation of both the computer simula- 
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Table 4,1 Variables recorded during simulation runs 


VARIABLE DEFINITION 


TIME 

DELA 

DELE 

DELR 

THRIL 

PA 

PADOT 

QA 

QADOT 

RA 

RADOT 

AXA 

AYA 

PSIA 

THEA- 

PHIA 

P 

Q 

R 

PDOTM 

QDOTM 

RDOTM 

AXCM 

AYCM 

PSIMB 

THEMB 

PHIMB 

XDDMB 

YDDMB 


time 

aileron deflection 
elevator deflection 
rudder deflection 
throttle input 
roll rate of airplane 
roll acceleration of airplane 
pitch rate of airplane 
pitch acceleration of airplane 
yaw rate of airplane 
yaw acceleration of airplane 
longitudinal acceleration of airplane 
lateral acceleration of airplane 
of airplane 
9 of airplane 
({) of airplane 

roll rate command to simulator 
pitch rate command to simulator 
yaw rate command to simulator 
roll acceleration measured on simulator 
pitch acceleration measured on simulator 
yaw acceleration measured on simulator 
longitudinal acceleration measured on simulator 
lateral acceleration measured on simulator 
^ of simulator 
0 of simulator 
<|) of simulator 

longitudinal acceleration of simulator without 
gravity component 

lateral acceleration of simulator without grav- 
ity component 
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tion of the motion and the actual simulator motion. This data 
was recorded at Langley on' their CDC 6600 computer. 

Figure 4,1 presents the aileron and rudder inputs to the 
simulation schemes, as previously shown in Chapter II. Table 
4.2 illustrates the four separate runs, and the data taken 
from each for use in the Ormsby model. Thus, there are twelve 
separate cases under evaluation. Both the rudder and the ailer- 
on inputs are simulated using the linear and nonlinear filters. 
For each of these four cases there are two simulated motion 
histories and one commanded motion history. 

The input to the Ormsby model is a subroutine known as 
STIM. The input to STIM is the time in seconds into the mo- 
tion history. This is computed in the main program. The out- 
put from STIM consists of three vectors - a specific force 
vector in g*s, a unit vector in the direction of gravity in 
g’s, and an angular velocity vector in radians/second. The - 
particular STIM subroutine used for this work can be found in 
the appendix. Basically, it reads the data from a file on 
disk in consecutive time order and places the desired data in 
the correct vector location. For example, when running the 
linear aileron roll data, the twentieth data item in the twen- 
ty-nine item list {see Table 4,1) is read into the first loca- 
tion of the angular velocity vector, after transforming it 
from an acceleration in degrees/second^ to a velocity in radi- 
ans/second, Thus, the STIM subroutine must be changed each 
time the model is run, to accomodate the new data. 
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\AXIS 
INPUT \ 

ROLL 

YAW 

AILERON 

-Linear 

Simulator 

(PDOTM) 

Nonlinear 

Simulator 

(PDOTM) 

Linear 

Simulator 

(RDOTM) 

Nonlinear 

Simulator 

(RDOTM) 

Command 

(PADOT) 

Command 

(RADOT) 

- 

RUDDER ■ 

Linear jNonlinear 

Simulator 1 Simulator 

! 

(PDOTM) . (PDOTM) 

; 


Nonlinear 

Simulator 

(RDOTM) 

i 

1 

r Command 

! 

■; (PADOT) 

< 

Command 

(RADOT) 


Simulator - recorded motions of the moving 
base simulator 

Command - requested motions of the moving 
base simulator made by the sim- 
ulation routine 


Table 4.2 Data used as input to Ormsby model 
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, The following four sections present the output of the 
model for the four major categories - aileron roll cues, ai- 
leron yaw cues, rudder roll cues and rudder yaw cues. 

4.2 Aileron Roll Cues 

Figures 4,2, 4,3 and 4,4 present the time histories of 
perceived angular velocity in response to aileron roll cues, 
using the linear and nonlinear washout schemes. In addition, 
the response to the commanded aileron roll is also shown. In 
each case, the perceived motion is approximately the same for 
the first thirteen seconds. The angular velocity rises grad- 
ually to a peak of .06 radians/second (3.5 degrees/second). 

This is consistent with the expected response to the 5*^/second 
input roll velocity of the pulse-type aileron cue. It is after 
this peak perceived velocity is reached that the interesting 
differences occur. 

But it is just at thirteen seconds when the second pulse 
is input. The linear and nonlinear washouts cause the perceiv- 
ed velocity to change direction, as indicated by the sign change. 
In the linear case, this change in direction does not occur un- 
til the end of the run, while in the nonlinear case it occurs 
at fifteen seconds. In both cases there is apparent confusion 
of direction. Just as there was in the first pulse, there 
should be a delay before the perceived angular velocity begins 
to return to zero. The experiment actually ends too soon, so 
the zero level is never reached. 




• Figure 4.2 Perceived angular velocity for simulated 


linear aileron roll cue input 




.0 
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A real difference can be seen in comparing the simulated 
cases with the commanded case. As can be seen in Figure 4,4, 
the commanded case behaves as predicted - there is a gradual 
increase to the maximum perceived angular velocity, and then a 
leveling off. Presumably, if the experiment had been carried 
past the second pulse, there would be a gradual return to zero 
in angular velocity 

In this case, then, the nonlinear filter acts to contain 
the confused perception involved in transferring the second 
pulse to the motion base. While it performs better than the 
linear filter, it presents motion cues which are not quite able 
to duplicate the desired motion perception. 
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4«3 Aileron Yaw Cues 

Figures 4.5, 4,6 and 4,7 present the perceived angular 
velocities output from the Ormsby model, for inputs of yaw 
cues for aileron motions. In .this case, the difference between 
the linear and the nonlinear washouts is evident. Again., the 
first thirteen seconds for each case are about the same - the 
expected response to a pulse input is the slow rise to a max- 
imum angular velocity, then a leveling off. This is the same 
response observed for the roll cues, as seen in Figures 4,2, 

4.3 and 4,-4. 

Thirteen seconds into the motion history, the second 
pulse is introduced. In the case of the roll cues, the 
motion transferred to the simulator was rather rough. But for 
the yaw cues, the simulation was very close to the desired mo- 
tion, This can also be seen by comparing Figure 2.9 with Figure 

2.10 - notice how smooth the nonlinear response is in Figure 

2.10 compared to the linear response in Figure 2,9, 

As before, the commanded motion to the simulator is smooth 
and presents, the expected response, A comparison of Figures 
4.5 and 4,6 shows that the nonlinear filter presented the sec- 
ond pulse with very little disturbance, while the linear filter 
caused a noticeable discontinuity in the motion. This is the 
anomalous rate cue which the pilots reported on in Table 2.1, 
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Figure 4.5 Perceived angular velocity for simulated 
linear aileron yaw cue input 
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4.4 Rudder Roll Cues 

Figures 4,8, 4.9 and 4.10 present the perceived angular 
velocities obtained from the rudder roll cue inputs. The sit- 
uation here is different from the previous aileron cases, sim- 
ply because the motion history in the rudder cue cases is much 
more complicated than in the aileron cue cases (see Figure 4.1) . 
It is not clear that the Ormsby model is equipped to handle 
such a rapidly varying motion history, and this must be kept 
in mind during an analysis. 

It does appear, however, that even in this more complex 
case, the nonlinear filter is able to contain the confused per- 
ceptions associated with transferring the pulse train to the 
motion base. Figure 4,10 shows that even the commanded input 
has wide motion discontinuity, which might lead to the conclu- 
sion that ^the Ormsby model has trouble handling this complex 
pulse train. Again, the perceived velocity gradually increases 
to a maximum, at about ten seconds. Had the experiment been 
continued past nineteen secohds , the' 'zero perceived velocity 
level would presumably gradually be reached. While there is 
some room- for argilment that the nonlinear filter better presents 
the motion cues in this case, it is a tenative argument at best. 
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4,5 Rudder Yaw Cues 

Figures 4,11, 4,12 and 4.13 present the final case - the 
perceived velocities obtained from rudder yaw cues. As in the 
previous case of rudder roll cues, the motion history is a com- 
plicated pulse-like train. But unlike the roll cues, the yaw 
cues seem to be transferred to the motion base more reliably, 
This was also true in the case of aileron inputs, - 

The motion histories for rudder ya;; cues are similar for 
the first ten seconds. This is attributed to the slow rise in 
angular velocity perception seen previously. The ten second 
rise time agrees with the rudder roll cue case. The nonlinear 
filter again does a better job of containing the discontinuous 
motion than does the linear filter. The commanded case is 
smoother than the simulated case, but the nonlinear filter 
does not change the commanded motion very much in the transfer 
to the motion base. 
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Figure 4,12 Perceived angular velocity for simulated 
nonlinear rudder yaw cue input 
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4 «6 Results 

The purpose of this investigation was to determine wheth- 
er or not there is a vestibular explanation for the results ob- 
tained by Parrish and Martin. These results (reported on in 
Chapter II) indicated that a nonlinear washout scheme provided 
better simulation fidelity than did the linear washout scheme. 
This result was not due to the fact that the nonlinear filter 
presented more of the motion cue; rather, it eliminated the 
false rate cue which arises in the use of the linear filter. 

In order to accomplish the goal of providing a vestibular 
explanation for the anomalous rate cue, the motion histories 
from the Parrish and Martin study were input to the Ormsby hu- 

r 

/ 

man dynamic orientation model. Included were aileron and rud- 
der motions with yaw and roll, cues, for each of the two wash- 
outs, The output from the model is the perceived angular vel- 
ocity of the pilot during the simulation. 

The outputs for each of the motion. -histories were present- 
ed in the preceding sections. Several results can be pointed 
out: 

• The yaw cues provide the most compelling case 
for a vestibular explanation. In the aileron 
yaw and the rudder yaw cases, the perceived 
angular velocities were "smoothed" considerably 
with the use of a nonlinear washout scheme as 
opposed to a linear washout scheme. The term 
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- "smooth? refers to the ability of the nonlin- 
ear filter to present a continuous motion clo- 
sely resembling the commanded motion, rather 
than the- 'discontinuous motion presented by the 
linear filter. The discontinuity which accom- 
panies the use of the linear filter has previ- 
ously been described as the fundamental differ- 
ence between the two filters - the anomalous 
rate cue. This false cue manifests itself in 

• the form of a jump in the perceived angular 
velocity of the pilot. 

• The results obtained for roll cue inputs were 
not so corroborative of the Parrish and Martin 
study as were the results for yaw cue inputs. 
They did, however, show some of the character- 
istics exhibited in the~^^w cue case. The non- 
linear filter contained the discontinuous jumps 
induced by the pulse train to a greater extent 
than the linear filter. The nonlinear filter 
was better able to transfer the commanded input 
to the motion base than the linear filter. 

This is evident in comparing Figures 4.11, 4.12 
and 4.13. 

• The explanation for the differences between the 
roll cues and -the yaw cues most likely could be 
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found in examining the mechanical differences 
between motion in the two axes; Intuitively, 
it can be argued that the yaw motion simulation 
(twisting about earth vertical) is an easier 
task mechanically than roll motion simulation 
(twisting about the horizontal axis) • No doubt 
a careful examination of the simulator base 
will reveal 'the cause of the differences obser- 
ved* 

• A comparison between the outputs for aileron 
and rudder inputs sheds some light on the use- 
fulness of the Oaonsby model. The aileron in- 
put consisted of two pulses, separated by thir- 
teen seconds, while the rudder input was a 
train of pulses. The Ormsby model has never 
been used with a complicated input such as the 
rudder input. But despite the fact that the 
output contains large motion discontinuities, 
even for the commanded case, it is still pos- 
sible to make a comparison between the linear 
and nonlinear schemes, and arrive at a conclu- 
sion similar to that reached in the aileron in- 
put case. Indeed, it does appear that the non- 


linear filter contains the discontinuous per- 
ceived angular velocity more effectively than 
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the linear filter. 

Thus it is seen that the Onasby model provides a vestibu- 
lar explanation for the sxibjectively acquired difference between 
the two washout schemes. The linear filter presents an anoma- 
lous rate cue as output from a pulse input, which the vestibu- 
lar system transforms into a discontinuous perceived angular 
velocity. The nonlinear filter does not present this false cue, 
and the resulting vestibular transformation provides a much 
"smoother" perceived angular velocity. In addition, the com- 
parison between the Orrasby model outputs from aileron and rud- 
der cue inputs gives insight to the model’s use as a simulator 
design tool. 
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CHAPTER V 
CONCLUSIONS 

This thesis began with a discussion of the importance of 
motion simulation in general, and went on to examine a parti- 
cular aspect of simulation - the washout filters used to con- 
strain the motion of the simulator and maintain the fidelity 
of the simulation. The two washout schemes examined here were 
a linear washout and a nonlinear washout. They differed in the 
types of filters used to washout translational cues. The lin- 
ear washout was seen to present a false rate cue in response 
to a pulse input, A subjective analysis of these two filters 
revealed that this false cue causes pilots to rate the fidelity 
of a simulation using the linear filter much lower than the 
same simulation using the nonlinear filter. 

Examination of physiological models of hviman dynamic orien 
tation led to the notion that such a model could be useful in 
comparing simulation schemes. The model used in this work, 
conceived by Ormsby, draws primarily on knowledge of the orien- 
tation infoanaation provided by processing information from the 
vestibular organs. Time histories for different motions were 



90 


input to the model in order to -evaluate the vestibular reac- 
tion to the linear ^d nonlinear filtering schemes* It was 
found that indeed the vestibular system reacts differently to 
the motion histories produced by the two filters. 

The next two sections present the conclusions of this work as 
they relate to the following two questions, first posed in the 
introduction : 

• Can the observed differences in simulation 
fidelity between the two filters be explain- 
ed using a physiological model of human dy- 
namic orientation? 

• What are the implications for this model as 
a drawing board tool in simulator design? 

The final section suggests avenues for further research in this 
area. 


5.1 The Vestibular Explanation Question 

Figures 5.1 and 5.2 present a recapitulation of figures 
shown in Chapter IV. They are the Ormsby model outputs for 
aileron yaw and rudder yaw cues, respectively, and they pre- 
sent the best cases for a vestibular explanation of the sub- 
jectively observed anomalous rate cues. In each case, the per- 
ceived angular velocity shows the expected gradual rise in re- 
action to the first acceleration in yaw. In the linear case, 
the second pulse (or pulses) causes discontinuities in the per- 
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Figure 5,1 Perceived angular velocity for three 
aileron yaw cue inputs. 
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ceived velocities. In the nonlinear cases, these disturbances 
are considerably reduced. This is best seen in comparing the 
simulated velocities with the commanded velocities. It can be 
seen that the nonlinear filter is better able to transfer the 
commanded motions than the linear filter. 

So, the physiological model enables a quantitative eval- 
uation of the differences in washout schemes to become a real- 
ity. It is now possible to 3cnow the outputs from the vestibu- 
lar sensors and to deduce a reason for the subjective ratings 
of the two methods. Indeed,, there is an anomalous rate cue 
sensed by the vestibular system - it manifests itself as a dis- 
continuous perceived angular velocity when the linear washout 
scheme is used, and that discontinuity is lessened considerably 
when the nonlinear scheme is used. 

The physiological model has performed the task demanded 
of it - it provided a vestibular explanation for the subject- 
ively observed differences between the two washout schemes. 

That difference was found in the differing perceived angular 
velocities which are the outputs from the model. 

While this was only a limited test of the perceptions in- 
volved in the motion simulation, it seems to validate the con- 
clusions reached in the Parrish and Martin study. It is also 
an additional validation of the model - since the predicted re- 
sponse to a pulse input is a gradual rise in perceived angular 
velocity to a maximum, and this is what was seen in every case. 
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the model appears to be functioning at a level consistent with 
available knowledge of the vestibular output. 

5,2 The Suitability as a Design Tool Question 

The question of the physiological model's appropriateness 
for use as a simulator design tool is a more difficult question 
to answer than the previous one. Certainly one could imagine 
the usefulness of such a model in simulation design. But the 
present case is a very limited one, and the small scope of this 
work should be taken into account in any conclusions which are 
drawn . 

Figures 5.3 and 5.4 present the roll cue.- input cases, as 
first presented in Chapter IV. The roll inputs did not propose 
as compelling a case for a vestibular explanation as the yaw 
inputs. But these figures are offered so that a comparison be- 
tween the aileron and rudder cases can be made. It is import- 
ant to remember that the inputs for the two cases are very dif- 
ferent - the aileron input is basically a pulse doublet, but 
the rudder input is a train of pulses. From this narrow inves- 
tigation it is hard to say whether the model really gives an 
accurate picture of the response to a complicated motion his- 
tory such as the rudder pulse train input. 

Assuming the model is proven to accurately portray the ves- 
tibular response to a complicated input, it appears that the 
model is applicable for simulation design purposes. In this 
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case, had the washout schemes been simulated on the computer 
rather than using a computer only to simulate the aircraft 
which in turn drives the motion base, the same motion histories 
could have been obtained. Then, without the necessity of set- • 
ting up an actual motion base, the same time histories could 
have been procured. Once input to the Ormsby model, the out- 
put would have shown the differences in simulation fidelity be- 
tween the two washout schemes. The same conclusions could have 
been reached ’ without ever having to use an actual mechanical 
simulator,' 

Thus, assuming the motion history of the part of the sim- 
ulator to be analyzed is sufficiently defined such that a com- 
puter simulation program can be written, the Ormsby model can 
predict pilot perceived angular velocities from that motion 
simulation. There is no need to use an, actual mechanical sim- 
ulator, and no; need to employ pilots for subjective analyses. 
The model is able to do the comparisons and predictions with 
confidence. 


0 '> 



98 


5.3 Suggestions for Further Research 

This work opens up several areas for further research: 

1 • It would be useful to understand how. certain 
parameters in each of the washout schemes affects the 
resulting motions of the simulator base, and the re- 
sulting perceived angular velocities of the pilot. By 
varying different parameters (such as in the preliminary 
filters or braking acceleration logic) new motion his- 
tories could be obtained. These, in turn, when input to 
the Ormsby model, could provide new insight into the 
workings of washout schemes. 

2. There are several revisions which suggest them- 
selves in regards to the Ormsby model. The necessity for 
tuning could be eliminated were the Kalman filters to be 
replaced by continuous Kalman filters, rather than the 
discrete filters currently in use. Also, more work should 
be done to verify that the model is indeed capable of 
handling complex motion histories. Finally, the model 
might be expanded to include visual and tactile cues, as 

..well as the vestibular cues it now employs, 

3. The model should be subjected to more rigorous 
tests of its ability to be used as a simulator design 
tool. One way which immediately suggests itself is to 
take a case such as the one examined here and do the 
testing in the opposite order. That is, run the motion 
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histories through the model, and then let the pilots do 
a sxibjective analysis. More extensive .use of the model 
will suggest areas for improvement, and begin to perfect 
it as a simulator design tool. 
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APPENDIX 

This appendix contains programming material used in the 
work presented in this thesis* The Man-Vehicle Leiboratory * s 
PDP 11/34 was the computer used for these FORTRAN programs. 

Most of the documentation for the main Ormsfay programs and as- 
sociated subroutines is taken from Borah [3] • 

A«1 Human Dynamic Orientation Model 

The listing which follows is the main module which imple- 
ments the Ormsby model of human dynamic orientation. Several 
changes have been made to the original program (the first three- 
by Borah) : 

1, Statements and routines which allowed for varying 
afferent base rates and additive random noise have 
been eliminated. Thus, all responses are average 
responses, and firing rates are those;. above the spon- 
taneous rate. 

2, statements were added to allow for non-zero long time 
constant, (variable TW in the program) , 

3, Comment cards were added for clarification. 
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4. Statements were added to calculate state transition 
matrices for any given update interval, for both 
canals and otoliths. Vectors TC, TPC, TO and TPO 
are no longer data entries, 

5. -.DATA statements replace data input cards. 

6. Kalman gains GKO and GKS were calculated for a .25 
second update interval, rather than the 1,0 second 
interval used by Oannsby and Borah. 

Table A.l lists the variables found in this program and 
their definitions. Several subroutines are needed to use this 
program and they are described in the next sections. Following 
the listing is a sample page of output. Table A, 2 describes 
the output variables seen on this page. 
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HLIHfttI OYK'f'niT.':. G!<IL"NTATION PROGR.^M 

n I HENS t DN RN P ( 1?0 L ) ^ TC ( 4 j 4 ) y TO < 3 y S ) ? TPC < 3 ? 3 ) y TPO ( 2 y 2 ) y DO C ( 3 ) y 
3 DU0(2> yCO<4> yCn<3) •/ THO ( 6 ) y XCH ( 4 ) . YCH ( 4 ) y ZCH ( 4 ) y CS < 3 ) » XOH < 3 ) y 

2 UluP (3) y t'OHso) yZr?H- y ( 3 ) y VO ( 3 ) y ZC < 3 > v XO ( 2 ) y YO C 2 ) - Or U ( 3 > ^ 

3 (3KS(3) yZ0(2) yCT0(3y3) »CT0(3y3) yUBRC<3) yOBROO) yA(3> yOC(3) y 

4 GP'C 4 ) y rWli < ) y rOH ( 3 ) y TBH ( 3 ) y TUiS ( 3 ) y T AS < 3 ) y AO ^ 3 ) y WO ( J ) y Y ( 3 ) y 

5 DOLB 1 3 ) y G: :f.U ( 3 ) y ElvO ( 3 > :• EUH ( 3 ) y EAH ( 3 ) y WPARE ( 3 ) y WPERP ( 3 W W7 i}T ( 3 ) y 

6 EAS ■: 3 ) y WSr- 0 ( 3 ) y WNCO ( 3 ) y WNCL C 3 ) y FN ( 3 ) v X < 3 ) y UO ( 3 ) y VG ( 3 ) y Vi- ( 3 ) 

DATA STATE! lEHTS 


DATA 

DATA 

DATA 

DATA 

DATA 

DATA 

DAT A 

DATA 
DA r A 
DATA 
DA1 a 
DATA 

data 

DATA 
DA I A 
IViTa 

data 
DA r A 
DATA 
DATh 
DATA 


DPR 757. 29370/ 

DT /.2500000000E 00/ 

NITP ;-B/ 

NDT 7240/ 

DVC /O V 224;l 30701----04 y 0 ♦ 47300;L20E-03 y 0 . 495801 SOE-027 
CC /-0.23578510E 02 y -0 . 1 1318880E 04 y -0 ♦ 637.L8550E 04 y 
0*A3;^,A19TOE 02/ 

ei.'t: /”0 -/T L91 7720E-03y 0. 0t41562OE-O5y 0 > 15035 L20E-03 y 
0*30559980E“01 ^ 

FSCC /1*570BOOOOE 00/ 


I'-AA 


00 / 

00 / 

/A . 47050^3oE“03y0»49053940E-02/ 


S'3C0 /'--O ♦ 7054000OE 00/ 


DVO 

on /I.SOOOOOOAE 03 . 0 . .OOOOOOOOE 04yO.OOOOOOOOE 00/ 
GKO /O * ci9/7505E-06 y 0 , 53 i 0665E--05 y 0 , 10602555E-02/ 
Gi-n /O- !.39'-i5i7E-05y04l0o2133E--04y0>21205tlE-02/ • 
FOTO /1,570B0E 00/ 

TC"0 /-’} ^3./-3E 00/ 

SO . u /- L , 57 0302 00/ 

SAOFAO /0-50000E 00/ 

OSPG /4.50000E 01/ 

DFAO /O.YAOOOi: 007 
XCii /A .OOwO-'.OOAE 00, 

0* 000000'.", OE 00/ 


O.O'-yOOOOOOE 00 y 0 » OAOVCOCOE OO^ 
DATA YOM /O.OOOOOOOOE 00 y 0 » OOOOOOOOE 00 y 0 , OOOOOOOOE 00 y 

O.OoOOOOOOE OOyO*OoOOOOOOE Oy, 


DATA 

data 

DATA 
DATA 
DATA 
DA TA 
DA ! A 
DA 1 A 
DATA 
DATi'i 
DATA 
riA 1 A 


0 >OOA0OOO'-''E 00''- 
ZCH /'0»OOOOOOOOE 00 :■ 
O.OOoOOOO'OE 00/ 
XOl-i /-O * Ob2h'OOOOE-0:l 
Y 0 i -i / 0 » 0 0 0 0 0 0 0 0 £ 0 ■ >) y 
ZDi ! /“O , 1 T325'400E"-01 


XC /O.OOOOOOOOE 
YC /0-0000<jOv'.'E 
ZC /0>0”00'’'000i- 
WNlVi •■'0 , OOO'"'',''' 

w {'( c L / u , On ‘ 0 ' > I ' 


XiJ /-O.lOGi 
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y - 0 , C3630000E-06 y -O ♦ 20000000F 00/ 
O»OOOO0dOOE 00y0,000e\.000i:. '.'0/ 
y -0 . 12650O00E-06 y -0 > 400000'.^OE 00/ 
00 y 0 , OOOOOOoOE 00 y 0 * 0000(}<n>O-i7 0 0 
00 y 0 ♦ OO00(.'OOOE 00 y 0 y OOOOOOOOE 00 / 

\ ■ 0 , ‘0 '0 0 0 O' 0 0 E O u y 0 ♦ k) 0 0 0 0 0 O 0 E 0; /' 

7-! 00 y 0 . OOnOOOC'OOE 00 y 0 * OOO'OOOOk.'F 'Oo/ 

^'0 7 0 . OOOOOO'OOE oo ■■ o , o'oooooooi;-. 00 / 

■ y-o. i332al :of:-oa/ 


Y U / 0 ♦ 0 « 0 0 0 0 < i •' 0 0 y 0 ♦ 0 < j 0 0 0 0 0 0 E 00/ 
on / -" ,2, 7.5.'..^90i£ -0i y O.ZAl'i !,460L--06/ 


Q-p poolt 
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C 

c 

c 


c 

c 

c 


c 

c 

c 


10 

15 


c 

c 

c 


20 


c 

c 

r. 


DATA DOLD /O^OOOOOOOOE 00 0 ♦ OOOOOOOOE OOy -1 ♦ OOOOOOOOE 00/ 
DATA UGFO /0*00000000£ 00 ? 0 » OOOOOOOOE 00 ? 0 ♦ OOOOOOOOE 00/ 
DATA UOF /O.OOOOOOOOE 00 r 0 ,OOOOOOOOE 00 >■ 0 ♦ OOOOOOOOE 00/ 
DATA TDUEL /3/oOOOOOOOE 01/ 

DATA TDPOS /6»OOOOOOOOE 01/ 

DATA TMC /0»25000000£ 00/ 

DATA FNOIBE /O.OOOOOOOOE 00/ 

SET UP FILE ASSIGNMENTS 

CALL ASSIGN(21? 'DK1U1BLT»N0NM 
CALL ASSIGN(22!- " DKl J RESUL f + SAR " ) 

WRITE(22!-5> 

FORMAT LINEAR AILERON ROLL 1 ') 


WRITE(22y700) DTyNITP 
WRITE(22y705) 


CANAL SPECIFICATIONS 

CALL STMCdiTyNITPyTPCyTC) 

DO 10 1=1 p4 

WRITE ( 22 y 710) TO ( I y 1 ) ? TC < I y 2 ) y TC ( I? 3 J ? TC (I y 4 ) 

DO 15 1=1 y 3 

WRITE ( 22 y 720 ) TPC ( I y 1) y TPC (I y 2 ) y TPC ( I y 3 X y D VC ( I ) 
WRITE(22y730) CC (1 ) y CC <2 ) y CC ( 3 ) yCCf4) 
WRITE<22y740) GKC ( 1 ) y GKC < 2 ) ^ GKC< 3 ) yGKC( 4 ) ' 

WRITE < 22 y 750 ) FSCC y TBCC y SSCC 
CALL EULER < FSCC y TSCC y SSCC y CTC ) 

WRITE < 22 y 765) 

OTOL I TH SPECIF I CATI OiNS 


CALL STMO (DTy NITP » TPO y TO ) 

DO 20 T=ty3 

WRITE<22y770) TO ( I r 1 ) y TO ( I y2 ) y TO ( I y 3) 


DO 25 1=1 y2 

WRITE ( 22 y 7B0 ) TPO ( I - 1 ) y TPO (I ■■ 2 ) y DVO C I ) 


WRITE(22y790) CO c I } yC0f2> yC0(3) 

WR I TE ( 22 y SOO ) GKO ( i ) y CK O ( 2 ) y OKO ( 3 ) y GKS (1 ) y (3KS < 2 ) y OKS { 3 > 
WR i; TE f 22 y S ! 0 ) F DTU y TOTO SOTO y S ACF i4C y OSI- U y DF nC 
CALL EULERCFO i 0 y fOTO y SO TO y C7 0 ) 

DO 27 1=1 y 3 


C 3 a)=CG( r^.KSACFAC 
WRlTE(22yS25) 


INITInLIlATIOij 


DO 32 r=i y4 

WRI TE (22 y 830 ) XCM ( I > y YCIJ ( I ) y ZCH < I ) 
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A(1)=XCH<4) 

A(2)=YCH(4) 

A(3)^ZCH(4) 

CALL C0TRN(A7CTC/1 ^WO) 

DO 36 1=1 y 3 

36 WRITE(22y836) XOH < I > y YOH (I > ? ZOH ( I ) 

A<1)=X0H(3) 

A(2)=Y0H(3> 

A(3)=Z0H<3) 

CALL COTRN ( A y CTO y 1 ? AO ) 
no 45 1=1 y 3 

45 WRITE<22yRS0) XC ( I ) y YC (I ) y ZC ( I) 

DO 50 1=1 y 2 

50 URITE(22yS55) XO (D y YO < I) y ZO (I ) 

WRITE<22y075) DOLDCl) pD 0LD(2) yD0LD(3) 

FP=1 » 0“EXP ( -DT/TDPOS ) 

FD=TDVEL 

WRITE < 22 y 880 ) TDOEL y TUPOS y TNC j FP y FNOISE 
. FNa)=EXP(-DT/TNC> ' 

FN ( 2 ) ( 1 ♦ ~FN ( 1 ) '> / DT-FN ( 1 ) 

FN(3)=1 ,~TNCH<(i*“'FN(l ) y/DT 

MAIN PROGRAM CYCLE 

FIND CURRENT STIMULUS IN HEAD COORDINATES 
(EVERY DT/NITP SEC»)J 

1* ANGULAR ROTATION VECTOR (TUH) AT (TIME)o 
2» SPECIFIC FORCE VECTOR <TAH) AT ( TIME4DT/2 ) , 
3. TRUE DOUN VECTOR AT (TIME+riT/2)» 

DO. 450 ITIME=lvHDT 
DO 100 I=lyNirP 
TIME= <ITIME-1 ) -KriT-f IH<DT/H1 TP 
CALL STIM(TIMEyTWHy TAHyTDH) 


TRANSFORM TO SENSOR COORDINATES 


CALL COTRN ( TWH y CTC y 0 y TUS ) 
CALL COTRN < TAH y CTO ? 0 y TAB ) 


SENSOR SriMULfiTION (EVERY DT/?LCTP 3EC*)i 

USING CURRENT STIMULUS VALUES y UPDATE STATE 
VECTORS FOR 3 CANALS (XCyYC AND ZC)y AfJD 
3 OTOI.ITHS (KOyVn AND ZO)y AND COMPUTE 
AFFERENT FIRING RATES ( CSX y CBY y CSZy OSX y OSY y OSZ ) » 


S=TWS(1 ) 


CALL SVUPD ( XC y TPC y IiVC s S y CSX y CC y 3 y 4 ) 
S=TWS(2) 

CALL SVUPD ( YC y TPC DVC y S - CS Y - CC y 3 y 4 ) 


S-=TWS(3' 

CALL SOU! 'll ZC y TT'C d VC y S y CS Z « CC y 3 y 4 ) 
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S=TAS(:!.) 

CALL SUUPD ( XO'y TPO ? DVO ^ S s- CSX 7 CO ? 2 y 3 ) 
S=TAS(2) 

CALL SOUPD < YO 7 TPO 7 DM0 7 S 7 OSY 7 CO 7 2 7 3 ) 
S=TAS<3) 

100 CALL SMUPD(Z07TP07riM07S70SZ7CS7273) 


OPTIMAL ESTIMATOR (UPDATE EMERY DT SEC*): 

GET CANAL AND OTOLITH SYSTEM STATE ESTIMATES FROM 
STEADY STATE KALMAN FILTERS* 

CALL SSKF < XCH 7 CSX 7 TC 7 CC 7 GKC 7 4 > 

CALL SSKF ( YCH 7 CS Y 7 TC 7 CC 7 GKC 7 4 ) 

CALL SSKF ( ZCH 7 CSZ 7 TC y CC 7 GKC 7 4 > 

CALL SSKFCXOHyOSXyTOyCOyGKOyS) 

CALL SSKF < YOH 7 OSY 7 TO 7 CO 7 GKO 7 3 ) 

CALL SSKF(Z0H70SZ7T07CS7GKS73) 

ENTER ROTATION RATE ESTIMATE MECTOR (CANAL ESTIMATE)* 


EWS(1)=XCH(4) 

EWS(2)=YCH(4) 

EWS(3)=ZCH(4) 

ENTER SPECIFIC FORCE ESTIMATE MECTOR (OTOLITH ESTIMATE)* 

EAS(1)=X0H(3) 

EAS(2)=Y0H(3) 


SACCULE NON-LINEARITY 

E AS ( 3 ) = AMAXl ( * 6 :;< ( ZOI-i ( 3 ) T * 4 1 69 ) - * 4 1 69 7 - * 4 ;i. 69 ) 

RESTORE MAGINITUDE OF OTOLITH ESTIMATE TO MALUE HELD 
BEFORE CONSIDERATION OF SACCULE NON-LINEARITY* 
(THEREFOREy NON-LINEARITY EFFECTS ONLY DIRECTION OF 
OTOLITH ESTIMATE) * 

CALL -NORMfEASyY) 

DO 130 1=173 

DUMMY=X0H ( 3 ) Jk^KZT-YOH ( 3 ) ^<>1<24Z0H ( 3 > 

EAS ( I ) =SQRT ( DUflHY ) ( I ) 

TRANSFORM TO HEAD COORDINATES 

CALL COTRN ( EWS 7 CTC 7 1 7 EWH ) 

CALL COTRN ( EAS 7 C TO 7 I ^ EAH ) 

PRINT STIMUl.US 7 SENSOR (^ND OPTIMAL ESTIMATOR MALUES* 
WRITE (22 7 900) TIME 

WRITE ( 22 7 91 0 ) TWl I ( 1) 7 TWS (1 ) 7 CSX 7 EUH ( 1 ) 7 TAH ( 1. ^ 7 TA5 ( 1) - 
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1 OSXyEAHCDjTriHd) 

WRITE 91.0)' TWH(2)dWB(2> 7 CSYi»EWH<2') ? TAH ( 2 ) ? tAS <2 ) ? 
1 OSYyEAHd) ?TDH<2) 

WRITE<22y910) TWH < 3 ) ? TWS ( 3 ) ? CSZ ? EWH (3 ) k TAH( 3 ) y TAS (3 ) >- 
1 0SZ7EAHC3) i-TEiH(3) 

WRITE (22^920) 

C 

C DOWN AND W ESTIMATOR < UPDATE EVERY DT SEC* ) ♦ 

C COMBINE OTOLITH AND CANAL ESTIMATES TO FORM ' 

C NEW ESTIMATE OF: 

C 1* PERCEIVED DOWN (DNEW) AT (TIME+DT/2)* 

C 2* PERCEIVED ACCELERATION (ACC) AT (TIME-fDT/2) * 

’C 3, PERCEIVED ANGULAR VELOCITY (WTOT.’' AT (TIME)* 

C 

CALL DOWN < BOLD ^ EWH ? EAH ? AO .» WSFO ? FD y DT ? TDPOS y DFAC y WOF y 
1 WNCOyWNCLyFN) 

450 CONTINUE 

C 

C FORMAT STATEMENTS 

C 


700 FORMAT ( // y " UPDATE INTERVAL:^' yF5 *2 y 'SECONDS* 'y 

1 'NUMBER ITERATIONS PER INTERVAL=dI3y//) 

705 FORMAT (//y' SEMI-CIRCULAR CANAL SYSTEM SPECIFICATIONS'?/ 

710 FORMAT< ' CANAL TRANSITION MATRIX=' y4E15*Q) 

720 FORMAT( ' CANAL SYS UPDATE MATRIX:^' y3E15*Sv 

1 ' CANAL DRIVING VECTOR^ ' y E15 * S ) 

730 - FORMAT(/y' CANAL SYS OUTPUT MATRIX^' ?4E15*8y/) 

740 FORMAT< ' CANAL SYS KALMAN GAINS ='y4E15*S) 

750 FORMAT</y CANAL ORIENTATION WRT HEAD PHI--= 'E12 *5? 

1 ' THETA=' -E12.5y ' PSI-=: ' y E12 * 5 ^ / ) 

765 FORMAT(////y ' OTOLITH SYSTEM SPECIFICATIONS '«/ ) 

770 >FORMAT( ' OTOLITH TRANSITION MATRIX==' y3E15 .8) 

780 FORMAT( ' OTOLITH SYS UPDATE MATRIX^' y2E15*3y 

1 ' OTOLITH DRIVING VECTOR”' y E15* 8) 

790 FORMAT(/y' OTOLITH SYS OUTPUT MATRIX=' y3£15*8) 

800 FDRMAT(/y' U.TR KAL GAINS=' y 3E12 ♦ 5y ' SAC KAL GAINS=' » 

1 3E12*5) 

810 FORMA-T(/y' OTOLITH ORIENTATION WRT HEAD PHI = '» 


1 


825 


FORMATC'l' 

830 


FORMAT ( ' 


1 

3E1 

836 


FORMAT ( ' 


1 

3E1 

850 


FORMAT ( ' 


1 

3E1 

855 


FORMAT ( ' 

875 


FOK‘MAT(/y ' 

880 


FORMA r ( / y ' 


E12*Sy' THETA--' yE12*5y ' PSI= ' y E12 *5 y / y 
" ■ SACFAC-='yE12*5v ' ■ 0 SENS PER G~' ? 

E12*Sy' 0 SYS GAIN ( DFAC ) - ' y E12 * 5 y / ) 

V/y ' SYSTEM 'INITIALIZATION' y/V) 

INITIAL STATE ESTIMATES* XCHy YCHyZCH=' y- 

. 8 ) 

INITIAL STATE ESTIMATES* XOH / YOH ^ ZGH"- ' •> 

. 8 ) 

TRUE CANAL STATE VECTORS XCyYCyZfd'y 
► 8 ) 

TRUF OTOLITH STATE VECTORS XO y YO v ZD'---- ' - 3F L5 * 8 ) 
DULD(ly2y3)=' y3E12*5y/) 

DOWN RA TE T"= ' y E12 * 5 y ' DOWf! POS T-'- ' y E12 > 5 y 


ORIGINAL PAGE IS 
OF POOR QUALITY 
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STOP 


900 

910 

920 


1 

2 

3 

1 


1 

n 


NON CONF T C0NS-'.Er2*5.-.P0s'ERRQR FAC~^ 


£ 12 . 5 ?/!- 
" 1 " ?//? ' 
format ( '■ 

'C 
'0 

FORMAT ( ' 
FORMATC' '■ 

'WPARE 


STOMAL MOISE FACTOR FiKTISE= 
SYSTEM SIMULATION' y/////) 


EST 


T='yF6*2y' W HB 
EST W HTi SF HB 
SF BQUN HB' 
' r9E12*5) 
RSCC , 
WPERP 

ACC' ) 


) 


U SENS 
SF SENS 


ROTO 

UTOT 


E12»5y/y 


C SIG 
0 SIG 


RPOS 
BNEW.' y 


STOP 

ENB 


/ 



108 


Table A*1 Variables used in main program 


AO (I) 

CC(I) 

cod) 

csd) 

CSX(I), CSY(I), CSZ{I) 

CTC(I,J) 

CTO(I,J) 

DFAC 

DOLD(I) 

DPR 

DT 

DVC(I) 

DVOd) 

EAH(I) 

EAS (I) 

EWH(I) 

EWSd) 

FD 

FN(I) 

FNOISE 

FOTO, SOTO/ TOTO 
FP 

FSCC, SSCC, TSCC 

GKC(I) 

GKOd) 

GKS(I) 


old otolith estimate 

canal sensor output 

utricle sensor output 

saccule sensor output 

current canal state vectors, sensor 
coordinates 

direction cosine matrix between head 
and canal, sensor coordinates 

direction cosine matrix between head 
and otolith, sensor coordinates 

steady-state gain of otolith estimate 

old DOWN value (=,46) 

degrees per radian (=57.29578) 

update interval for DOWN . estimator 

canal sensor driving vector 

otolith sensor driving vector 

current otolith specific force esti- 
mate, head coordinates 

current otolith specific force esti- 
mate, sensor coordinates 

current canal angular velocity esti- 
mate, head coordinates 

current canal angular velocity esti- 
mate, sensor coordinates 

« TDVEL 

constants for first-order, filter 

signal-to-noise factor 

Euler angles for head and otolith sen- 
sor coordinate transform 

position error factor 

Euler angles for head and canal sen- 
sor coordinate transform 

canal Kalman gains 

otolith (utricle) Kalman gains 

otolith (saccule) Kalman gains 


ORIGINAL MGE is 
OF . POOR QUALM 
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Table A.l 

NDT 

NITP 

OSFG 

0SX{1) , OSY(I) , 

SACFAC 

TAH(I) 

TAS(I) 

TC(I,J) 

TDH(I) 

TDPOS 

TDVEL 

TIME 

TNG 

TO(I,J) 

TPC(I,J) 

TPO(I,J) 

TWH(I) 

TWS(I) 

WNCL(I) 

WNCO(I) 

WO (I) 

WOF(I) 


continued 


length of motion history 

number of sensor updates per DT 

otolith afferent firing rate per g 

OSZ(I) current otolith state vectors, sensor 
coordinates 

saccule factor (=.5) 

current stimulus specific force, 
head coordinates 

' current stimulus specific force, 
sensor coordinates 

state transition matrix for canal 
Kalman filters 

DOWN 

60 second time constant for DOWN 
position 

35 second time constant for DOWN 
angular velocity 

current time in seconds 

,25 second time constant for uncon- 
firmed canal estimate 

state transition matrix for otolith 
Kalman filters 

state transition matrix for canal 
sensor update 

state transition matrix for otolith 
sensor update 

stimulus angular velocity, head 
coordinates 

stimulus angular velocity, sensor 
coordinates 

low frequency portion of WNCO 

previous unconfirmed canal angular 
velocity estimate 

old otolith estimate 

low frequency portion of otolith 
angular velocity estimate 



Table A.l 


WSFO(I) 

XC(I) , YC(I) , ZC(I) 
XCH(I), YCH{I), ZCH(I) 
XO(I) , YO(I) , ZO{I) 
XOH(I), YOH(I), ZOH(I) 
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previous otolith angular velocity 
estimate 

old canal state vectors, sensor 
coordinates 

current canal state vectors, sensor 
coordinates 

old otolith state vectors, sensor 
coordinates 

current otolith state vectors, sensor 
coordinates 
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Rsc:c ROTO RPOR UPARE 
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T= 0.50 W HD U SENS C SIG C EST U HD 

-0. 13216E-02-0 , B0697E-03-0. 2665SE-01-0 . 1 1467E-02 
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Table A. 2 Variables output from model 


W HD 
W SENS 

C SIG 

C EST W HD 

SP HD 
SF SENS 
0 SIG 

0 EST ST 

DOWN HD 

RSCC 

ROTO 

RPOS 

WPARE 

WPERP 

WTOT 

DNEW 

ACC 


angular velocity vector, head coordinates 

angular velocity vector, sensor coordi- 
nates 

canal signal; afferent firing rate from 
three canals 

canal estimate of angular velocity vector, 
head coordinates 

specific force vector, head coordinates 

specific force vector, sensor coordinates 

otolith signal; afferent firing rate from 
three otoliths 

otolith estimate of specific force vector, 
head coordinates 

unit vector in direction of. gravity, head 
coordinates 

canal contribution to DOWN, head coordi- 
nates 

otolith contribution to DOWN, head coor- 
dinates 

rotation vector to null difference between 
SF and DOWN 

angular velocity perception parallel to DOWN, 
head coordinates 

angular velocity perception perpendicular to 
DOWN, head coordinates 

total perceived angular velocity, head coor- 
dinates 

perceived DOWN vector, head coordinates 

perceived acceleration vector, head coordi- 
nates 
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A, 2 STibrou-tine STIM 

Subroutine STIM is the 'stimulus' routine called by the 
main program* It is this program which is altered for differ- 
ent applications of the model. The particular subroutine 
listed here is the one used in this thesis research. Basic- 
ally, it reads the data from the Langley motion histories from 
a file on a disk. The desired angular velocities are placed 
in the proper angular velocity vector locations by this pro- 
gram, Table A, 3 lists the filenames, data locations and vector 
locations for each of the twelve cases exetiained. Table A, 4 
lists the variables used in this program, and their definitions. 

Note that the STIM routine must return staggered angular 
velocity and specific force values, as required by the main 
program. The value of W must correspond to time T, while the 
values of A and D must correspond to time T+DT/2. This require- 
ment is illustrated in Table A. 5, which also gives the print- 
out times for the variables which are output from the main pro- 


gram 
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10 

SUBROUTINE STIH(Ti>Uj.AyD> 
niMEWSrON WC5) vA(3) ?D(3>rIi8TA<4J 
DFR-57 ,2957B 
BO 10 1=1 ?o 

wai=o*o 
D < I ) =0 1 0 
A(I)=0.0 


A<3)=”1»0 
D<3)=-1 *0 
DO 20 1=1? 5 

READ<2a.?15?ERR=i6> (IiATA< J> ? J=1 ?65 

15 

F0RfiAT(lX?6El-3.6> 


GO TO 18 

16 

WRrTE(7?15> CDATA(K)7K»1?6> 

18 

1F(X*£Q*2> W Cl ) =DA TA < 1 ) iKT/BPR 

20 

CONTINUE 

STOP — 

RETURN 

END 


OEIGINAL VA.GE IS 
OS' PCiOR QiOAiaTY 
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Table A. 3 Twelve test cases used by STIM 


FILE NAME 
DBLT.LIN 

DBLT.NON 

DBLT.LIN 

DBLT.LIN 

.DBLT.NON 

DBLT.LIN 

RDDDR.LIN 

RUDDR.NON 

RUDDR.LIN 

RUDDR.LIN 

RUDDR.NON 

RUDDR.LIN 

yaw 


DATA NAME & 
LOCATION * 

ANGULAR VELOCITY 
VECTOR LOCATION 

PDOTM 

(4,2) 

1 

PDOTM 

(4,2) 

1 

PADOT 

(2,1) 

1 

RDOTM 

(4,4) 

3 

RDOTM 

(4,4) 

3 

RADOT 

(2,5) 

3 

PDOTM 

(4,2) 

1 

PDOTM 

(4,2) 

1 

PADOT 

(2,1) 

1 

RDOTM 

(4,4) 

3 

RDOTM 

(4,4) 

3 

RADOT 

(2,5) 

3 


CASE 

Simulated aileron 
linear roll 

Simulated aileron 
nonlinear roll 

Commanded aileron 
roll 

Simulated aileron 
linear yaw 

Simulated aileron 
nonlinear yaw 

Commanded aileron 
yaw 

Simulated rudder 
linear roll 

Simulated rudder 
nonlinear roll 

Commanded rudder 
roll 

Simulated rudder 
linear yaw 

Simulated rudder 
nonlinear yaw 

Commanded rudder 


* ' Data ‘ location taken from a 5 X 6 macrix of variables 
listed in' Table 4,1 
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A(I) 

D(I) 

T 

W(I) 


Table A. 4 Variables used in STIM 

stimulus specific force, head coordinates 

unit vector aligned with gravity, head coor- 
dinates 

current time in seconds 

stimulus angular velocity, head coordinates 
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Table A, 5 STIM variables and printout variables [3] 


VARIABLE 

COORDINATE 

FRAME 

UNITS 

STIM - 
VALUES 

T+DT/2 

COMPUTE 
AT TIMES 

T 

PRINTOUT 
VALUES AT 

T+DT/2 

TIMES 

T 

A 

head 

g 

X 




D 

head 

g 

X 




W 

head 

rad/sec 


X 



W HD 

head 

rad/sec 




X 

W SENS 

sensor 

rad/sec 




X 

C SIG 

sensor 

ips 




X 

: EST W 

head 

rad/sec 




X 

SF HD 

head 

g 



X 


SF SENS 

sensor 

g 



X 


0 SIG 

sensor 

ips 



X 


0 EST SF 

head 

g 



X 


DOWN HD 

head 

g 



X 


RSCC 

head 

rad 




X 

ROTO 

head 

rad 




X 

RPOS 

head 

rad 




X 

WPARE 

head 

rad/sec 




X 

WPERP 

head 

rad/sec 




X 

WTOT 

head 

rad/sec 




X 

DNEW 

head 

g 



X 


ACC 

head 

g 



X 
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A»3 Subroutine DOWN 

. Subroutine DOWN implements the logic for determining the 
perceived direction of gravity and the perceived angular vel- 
ocity. Figures 3,7 and 3.8 illustrated this logic, and it was 
discussed in Chapter III. Table A. 6 provides the list of var-r 
iables used in the subroutine along with their definitions. 
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C 

c 

c 

C- 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 

c 


10 


11 

12 


A **y 


SUBROUTINE DOOM ( BOLD i- UN SN ? SO ^ WSFO » T yDT »TDP3 y 
1 BFACjWOF?WNCO?UNL,FN) 

DOWN ESTIMATOR AND U ESTIMATOR 

DOWN IS DETERMINED BY RELYING ON LOW FREQUENCY 
OTOLITH ESTIMATES y CANAL ESTIMATES WHICH ARE 
CONSISTENT WITH HIGH FREQUENCY OTOLITH ESTIMATES 
AND HIGH FREQUENCY PORTION OF CANAL ESTIMATES 
NOT CONFIRMED BY OTOLITHS* 

W IS DETERMINED BY CANAL ESTIMATES PARALLEL TO 
DOWNS ROTATION RATE OF DOWN? AND HIGH FREQUENCY 
PORTION OF CANAL ESTIMATES PERPENDICULAR TO 
DOWN MINUS ROTATION RATE OF DOWN* 

DIMENSION BOLD <3 ) ?WN(3> ?SN<3) ? SO < 3 ) ? WSFO ( 3) ? F ( 3 ) ? 

1 WSF<3) ?X(3) ?W0F(3) ?RSCC<3) ?R0T0(3) ?RT0T(3) ? 

2 DNEW(3) ?RP0S<3) ?WPERP<3) yDA0G<3> yWPARE(3) ? 

3 WT0T<3) ?ACC<3) rANG(3) vW0D<3)?W0DN(3) ? 

4 WNCO (3 > ? WNL < 3 ) ? FN ( 3 ) ? WNC ( 3 ) ? WNCH ( 3 ) ? HROTO ( 3 


SFMAG^SQRT ( SN (1 ) :KJK2+SN ( 2 ) >K^:2+SN ( 3 ) -Y>K2 ) 
FP0S=l*0-EXP<-< (SFMAG/DFAC))K>K( *25) ):KDT/TDPS) 
TDVEL=:T 

F(1.)=EXP(-DT/TDU£L) 

F ( 2 ) ==TDyEL5K (1 * -F ( 1 ) ) /DT-F ( 1 ) 
F(3)=1*-TDVEL>:?<1*-F<1) )/DT 
CALL CROSS(SOsSN?WSF) 

CALL NORMCWSFsX) 

CAL i . UANG ( SO ? SN ? AL'GSF ) 

DO 10 1=1 7 3 
WSF(I)=ANGSF>!<Xa) 

WQF ( I ) =F < 1 ) 5KW0F (I ) -PF ( 2 ) ^WSFO ( I ) +F ( 3 ) :!<WSF < I ) 
WOD(I)=WSF(I)-'WOF(I> 

WSFOa)=-WSF(I) 

WODM=SORT ( WOD ( 1 ) j;-j;K2+WGD ( 2 ) :{c:!<2+W0i.i ( 3 ) tr2 ) 

CALL NORM ( WOD ? WODN ^ 

WCPWD=WH ( i ) :!<WnON ( 1 > -P WN < 2 ) WODN 1 2 ) PWN ( 3 ) ::<WODN ( 3 ) 

IF(WCPWD) 12? 12? 11 

WCPWD=0 * 0 

WMAiG=~WCPWDY'0 T 

no 13 1=173 

X ( I ) =WODN f [ ) ;KANIN:I ( WHAG 7 WODil > 

WNC 1 1 ) =-WN ( I ) Ml\ r-x < I ) 

WNL ( I ) "FN < 1 ) ^:WNL < I ) -PFN ( 2 ) TWNCO ( I ) +FN ( 3 ) .KWNC (I)- 
WNC0<1)=WNC( !■> 

WNCH ( i ) --=WNC ( I ) -WNL ( I ) 

X( J>=X(I)iWNCl-KI'i 
C A L L C H 0 G S Ui 0 L D ? S N ? 0 T 0 ) 

CALL NORM-: ROTO ?F> 
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STO? 


15 


20 


30 


40 


45 

50 

100 


UiOP,‘iU’it-=WOi" .i. ) :n~ < .1. ) 4W0F ( 2 ) :KF ( 2 ) -f WOF ( 3 ) ( 3> 

BD 15 


ROTO ( J ) ==VvOPf-, KH.'kF ■, I ) 
HRQT0(l>==i?0)0<i:)/2» 

CALL ROTA i £ < UULD y H'‘''nTn y F ) 


WCPAFci-i- 


' ( '! 


< U +A < 2 ) ;:<F ( 2 ) +X ( 3 ) tF ( 3 ) 


DO 20 T~l. "3 

RSCC (J ) =: C X !. r ) -Wt.'Pf-'sRM-M- ( I) ) 

RTOT ( X J "RSCC < I ) -i-RO 10(0 
CALL ROTATF < BOLD ? R~\ (3T ^ DNEW ) 

CALL K'iVAB ( DNCW s- SN v FEE ) 

PHX=FPOS-|tFI£E 

CALL CROOG!;XC:r:i;?SN!>RPnSJ 

CALL NOF*;H(RPas?X) 

DO 30 T~1r3 

Rpos<i)==rHn:<xa > 

CALL ROTATE < HNEl'l ^ RPOS ? X ) 

CALL' HOF^fHXyriHEl'!) 

DO 40 1=1? 3 
x(i)=DnLo<n hBMEwa) 

CALL NurCi(X?DiiO0? 

WPARM=WtX C 1 ) 4<DAVG ( 1 )+WN ( 2 ) ^DAOG ( 2 ) +WHX 3 ) :{cDAVG (3 ) 
CAL L CROSS ( DHEIJ BOLD ? WPERP ) 

CALL yPERP?X) 

CAi L •Aii-iGi::!nLri?riNEw?PH n 
DO 50 .i; = J.?3 


UPERP ( X ) -::X ■. T :s MRH I /DT 
WPARE ( X ) ==OPARh>KDAVG < X > 
DOLD(T)-^Ti;:OK : ^ 
S0(i)=SF!n;> 

UTOT < I ) - iC: FIR? ( X y -fMEARF 


ACC(I)=Dr 

yRXTE(22? 


ACTiiPFl-J<; I. )“bNCI) 

X 00< ^ RGGC >: X > V ROTO (I ) ? RPOS ( I > 
WTOT (1 ^ ? DNEl<{ ( I ) ? ACC ( I ) 


» yPARE a ? 


CONTINUE 


FORhA'K-' ^?SE12.5) 

RETUF5M 

END 


WPERP ( I ) ? 
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Table A. 6 Variables used in subroutine DOWN 


ACC 

AN6SF 

DAVG(I) 

DFAC 

DNEW(I) 

BOLD (I) 

DT 

FEE 

FN(I) 

FPOS 

HROTO(l) 

PHI 

ROTO (I) 
RPOS(I) 

RSCC(I) 

RTOT(I) 

SPMAG 

SN('I) 

sod) 

T 

TDPS 

TDVEL 

WCPARM 

WCPWD(I) 

WMAG 

WN(I) 

WNC(I) 

WNCH(I) 

WNCO(I) 


perceived acceleration vector 
angle between SO and SN 

unit vector in the direction of BOLD + DNEW 
steady-state gain of otolith estimate 
current value of DOWN 
old value of DOWN 

update interval for DOt-JN estimator 

angle between DNEW. and SN 

constants for first-order filter 
^-DT/TDPS 

ROTO/2 

FEE*FPOS 

component of WOF perpendicular to SN and BOLD 

rotation vector to eliminate integration 
errors 

canal contribution to DO^TN 
RSCC+ROTO 

specific force magnitude 
current otolith specific force estimate 
old otolith specific force estimate 
=TDVEL 

time constant for DOWN position 

time constant for DOWN angular velocity 

magnitude of canal angular velocity estimate 

canal estimate parallel to high frequency 
angular velocity otolith estimate 

angular velocity magnitude 

current canal estimate of angular velocity 

current canal estimate of angular velocity 
not confirmed by otolith estimate 

high frequency portion of WNC 

old WNC 
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Table A. 6 

concluded 


WNL(I) 

low frequency portion of WNC 


WOD 

WSF - WOF 


WODM 

magnitude of WOD 


WOF(I) 

low frequency portion of WSF 


WOPARM 

magnitude of otolith angular velocity 

estimate 

WPARE(I) 

canal angular velocity parallel to DNEW - 

WPARM(I) 

canal angular velocity perpendicular to DNEW 

WPERP(I) 

system angular velocity perpendicular 
DNEW 

to 

WSF{I) 

angular velocity of otolith estimate 


WSFO(I) 

old WSF 


WTOT(I) 

WPARM + WPERP 
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A»4 Subroutine library 

The next listing contains the rest of the subroutines 
used by the model. They are mostly self-explanatory. Sub- 
routines STMO and STMC are new routines designed to calculate 
state transition matrices for any given update interval. They 
implement the equations for the systems described in Tables 
A. 7 and A. 8. 
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SUBROUTINE SUUPD (Xs-T?D!-Si-YrCi-N!>M) 

C 

• C STATE VECTOR UPDATE': 

C 

C X(NEU)=T>i'X<QUD + D;SS 

C Y(NEW)=C>}{X<NEW) 

C 

C WHERE 

C N IS DIMENSION OF STATE UECTDR 

C X IS STATE VECTOR 

C T IS TRANSITION MATRIX 

C n IS DRIVING VECTOR 

C S IS STIMULUS 

C Y IS OUTPUT (AFFERENT FIRING RATE) 

C C IS OUTPUT MATRIX 

C 

DIMENSION X(N) y rCNi-N) jD(H)?C(M) ?R<9) 
DO 5 I=5l,N 
5 R(I)=Xa) 

DO iO 1=1 yN 
X(i)=D(I)>i<S 
DO 10 J=lyN 

10 X(I)=Xa)-fT(IyJ);i<R(J)' 

Y=C(M)JitS 
DO 20 1=1 yN 
20 Y=Y+C (I ) $X (I ) 

RETURN 

END 


SUBROUT I NE NORM ( A y AN ) 

C 

C AN = UNIT VECTOR IN DIRECTION OF VECTOR A 

C 

DIMENSION A(3)yAMC3) 

AM=SQRT < A ( 1 ) jK;!<2+A ( 2 ) in2-\-P\ ( 3 ) ) 

IF(AM-1 *E-06) 5y5y<S 

5 AM=t*E-OA 

6 BO 10 1=1 y 3 

10 ANC I)=A(I)/AM 

RETURN 

END 


C 

C 

C 

C 

C 

C 

c 


SUBROUTINE SSKF ( Xl-I y Y y TM y C y GK y N) 

STEADY-STATE KALMAN FILTER (UPDATE EVERY DT SECONDS) 
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C TM IS TRANSITION MATRIX 

C GK IS KALMAN GAIN MATRIX 

C Y IS SENSOR SYSTEM OUTPUT ' 

•C C IS OUTPUT MATRIX 

C 

DIMENSION Xi-l^N) y TM<NfN) yC<N) i-GK<N) :>S(9) 
DO 40 I=iyN 
SCI):=0»0 
DO 40 J=lyN 

40 S(I)=S(I)+TM(I y J) 

EM=0*0 
DO 45 1=1 yN 
45 EM=EM+S<T)>!cCa) 

DO 50 J=lyN 

50 XH( J)=S< J)+GK( J)H<(Y~EM) 

RETURN 

END 


SUBROUTINE ROTATE ( A? R y AR ) 

C 

C .AR = A ROTATED ABOUT R BY AN ANGLE (RAD) 

C EQUAL TO THE MAGNITUDE OF R 

C 

DIMENSION A(3) yR(3) yAR(3) yAP(3) yAPN(3) 
CALL CROSS(RyAyAP) 

CALL NORM(nPyAPN) 

AMAG=SQRT ( A (1 ) i!<:K2+A < 2 > 55C5!<2 f A ( 3 ) ^^^2 > 

PH I=SQRT ( R (1 ) >K5K2+R ( 2 ) 2+R ( 3 ) tt2 ) 

DO 10 1=1 y 3 

10 • AR<I)=AMAGJ{(SIN(PHI>;S{APN(I )+COS(PHI)>!<A(I) 

RETURN 
END 


SUBROUTINE COTRN ( A y B y N y C ) 

C’ 

C COORDINATE TRANSFORMS 

C FROM HEAD TO SENSOR Jl- N=0 

C FROM SENSOR TO HEAD IF N--=I 

C 

C A = ORIGINAL VECTOR 

C C = TRANSFORMED VFCTOR 

C E = transformation MATRIX 

C 
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fUMENSION At3) ?E(3?3) ?C(3> 

IF<N) i0yl0r20 
10 no 15 1=1 v3 

15 C < I ) =B < I !- 1 ) A ( 1 ) +B (I ? 2 ) < 2 ) 4-B < I J-.3 ) >icA ( 3 ) 

GO TO 30 

20 no 25 1=1 ”3 

25 C(I)=B(1 7I)>KA(1)+B(2jI)>KA(2)+B(3»I)>KA<3) 

30 RETURN 

END 


SUBROUTINE OANG ( A ? D ? PHI ) 

PHI = ANGLE BETWEEN A AND B 

DIMENSION A(3) yB<3) j AN(3) ?BN(3) 

CALL NORM ( A AN.' 

CALL NORM <3? BN) 

X=AN ( 1 ) )!<BN ( 1 ) -f AN < 2 ) Jl<BN ( 2 ) +AN < 3 ) ^<BN ( 3 )■ 
IF(X.GT*1*0) X=1*0 
Y=SQRT(1 »-X«>K2) 

PHI=ATAN2<Y?X) 

RETURN 

END 


SUBRul IT I NE Ci’^OSS < A ? 3 y C ) 


C = A X B 


D I HENS I ON A ( 3 ) ? B ^ 3 ) " C < 3 > 
C ( 1 ) = A < 2 ) 5:. B < 3 ) - A ( 3 > >|i B ( 2 ) 
C(2)=A(3)':<S(i)-ACl)5:sE(3) 
C ( 3 ) =A ( 1 > im < 2 ) C 2 ) < 1 ) 

RE TURN 
END 


SUBROUTINE EUI..ER ( F * T y S ? CT ) 
r 


OF POOR 


PAGE IS 
QUALtTY 
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C PRODUCE DIRECTIGH COSINE MATRIX <CT) GlO-EN 

C EULER ANGLES (F?T AMD S) . 

C 

DIMENSION CT<3r3) 

C7<lyl)=C0S<S)«C0S(F>-C0S<T)>KSIN<F)Jf(Si^N<S> 
CT<2yl)=-<SIN<3)4tC0S(F)+CaS<T)#SIN<F).'{cC0S'<S> j 
CT(3fl)=SIN<T)»SIN<F> 

•- CT<l72)=C0SCS)>S(SIN(F)+C0S(T)*C0S<F)«SIN<S> 
CT<272)=C03<T):^{COS<F)J|{COS(S)-SIN(S):{«SIN<F) 
CT<372)=-SIN<T)#C0S<F) 

CT<l73)=SIN(T)JtcSIN<S) 

CT<273)=SIN(T)1KC0S<S5 
. CT(373):=caS(T> 

.DO 10 I=lf3 

10 URITE<22»100) CTdf 1) >CT(l72) »CT<1»3) 

100 FQRNAT( ' CT= ' 7 3E15 * 8 ) 

RETURN 

END- 


SUBROUTINE STM0<DT7NITP7TP0»T0) 

C 

C SUBROUTINE TO CALCULATE OTOLITH STATE 

C TRANSITION MATRICES 

C DT - UPDATE INTERVAL FOR KALMAN FILTER 

C NITP ~ NUMBER OF SENSOR UPDATES PER DT 

C TPO - 2X2 STM FOR OTOLITH SENSOR UPDATE 

e TO - 3X3 STM FOR OTOLITH KALMAN FILTER SYSTEM 

C 

DIMENSION TP0<2r2) 7T0<373) 

C CALCULATE TPO 

T=DT/NITP 

TP0<l7l) = ^200♦5^cEXP<-♦2^}{T)-♦2:fcEXP<-200♦^^cT) >/199,S 
TF0<27i) = (EXP(-*2:^T)~EXP(-200*--KT) )/19?»S ‘ 

TPO (17 2) ( EX? < -200 . ^T ) -EXP < - * 2^7 ) ) /199 . 8 

TPO(272) = <200*-lcEXPC-200.:icT)-*2»EXP<-.2SiT) )/199*3 
C CALCULATE TO 

T=DT 

TO < 1 7 1 ) = < 200 ♦ :tcEXP ( - . 25i«T > - ♦ 2JICEXP < -200 . 2i:T ) ) / 1 99 ♦ S 
T0<27l)=(£XPC-.2YT)-EXP(-200*>KT) )/199.3 
T0(37l) = (EXP(-200»^'«T) )/39760»2+(EXP(-«2.=itT) )/159«S4 
1 -iEXP<-T) >/159.2 

TO ( 1 7 2 ) --^ 0 ^ ( EXP ( -200 * ^KT ) -EXP ( - ♦ 2:KT ) ) /1 99 » 8 
TO<272> = (200..«EXP<-200,-'fcT)-,2:<cEXP<-»2ST) )/199*3 
T0(372) = ^.EXF(~T) )/lS9*2-( ♦ 2:icEXP ( - . 2^T ) ) /159 * 34 
1 - < 200 » .'f^EXP < -200 ♦ :5« T > ) /39760 ♦ 2 

T0<l73)=0*0 
'70(273)^0.0 
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TQ(3»3)=EXP<-n 

RETURN 

END 


SUBROUTINE STMC ( DT y NITP y TPCy TC ) 


SUBROUTINE TO CALCULATE SEMICIRCULAR CANAL 
STATE TRANSITION MATRICES 

DT - UPDATE INTERVAL FOR KALMAN FILTER 
NITP - NUMBER OF SENSOR UPDATES PER DT 
TPC ~ 3X3 STM FOR CANAL SENSOR UPDATE 
TC - 4X4 STM FOR CANAL KALMAN FILTER SYSTEM 


C 


C 


I 

1 

1 

1 


1 

0 

A. 

1 
O 


I 


DIMENSION TCC4y4) rTPC(3y3) 

CALCULATE TPC 
T=DT/NITP 

TPC ( 1 y 1 ) = < ♦ 0557SJicEXP ( - ♦ 03322:KT ) - ♦ 03322>KEXP < - ♦ 055765ST ) ) / ♦ 02254 
TPC ( 2y 1) =-200 ♦ 033.YEXP ( - » 0S576:KT ) /^ ♦ 50674 

+200 * 0555S+EXP ( - ♦ 03322>!<T ) /4 * 50725 
+ ♦ 089>f!EXP < "199 + 9998>KT ) 739982 ♦ 118 
TPC ( 3 y 1 ) =EXP (-199. 9998;!(T ) 739982 . 1 1 8+EXP < - ♦ 03322-KT ) 74 . 50725 
-EXP (-♦ 05576>i^T) 74 .30674 
TPC (1 y 2 ) =- . 370373ic ( EXP (-199. 9998:|<T ) 739982 * 1 1 8+ 

EXP ( - . 03322JKT ) 74 . 50725-EXP C - . 05576>i^T ) 74 . 50674 ) 

TPC ( 2 y 2 ) =2 . 47492?KEXP < - . 05576«T ) - 1 . 4744S>KE:XP ( - ♦ 03322;KT 5 
- . 0004452^EXP (-199. 999SAT ) 


TPC ( 3 y 2 ) = . 05576:KEXP ( - . 05576»T ) 74 . 50674 

- ♦ 03322*EXP ( - . 03322>!<T ) 74 . 50725 

-199 . 9998:{CEXP( -199 > 9998.’KT > 739982 . 118 
TPC ( 1 y 3 ) = . 37037:-}: ( . 03322^«EXP ( - . 03322>!<T ) 74 . 50725 
+199 . 9998+XrXP (-199 . 999SYT ) 739982 .118 

- * 05576H4EXP ( - . 05576>!<T ) 74 ♦ 50674 ) 


TPC ( 2 y 3 ) =-l 7 . 7966'^ ( ♦ 03496JKEXP ( - . 05576>KT ) 74 . 50674 
•- . 01242>!cEXP ( - . 03322JKT ) 74 . 50725 


-1 99 . 979$EXP (-1 99 .9998*T ) 739982 . 1 18) 


TPC ( 3 y 3) =1 . 000445.'KEXP ( -199 . 9998>HT > + . 00024484JXEXP ( 


03322.t r ) 


.05576>KT) ) 7. 02254 


- . 0006899>!<EXP ( - . 05576iKT ) 

CALCULATE TC 
r=DT 

TC (1 y I ) = ( .05576:{{EXP (~ . 03322>KT ) - >03322^!iEXP ( • 

TC ( 2 y 1 ) =-200 . 033^KEXP ( - . 0 5 5 765K I’) 74 . 50674 

+200 . 05358;!iEXP ( - . 0332?;i<‘l ) 74 ♦ 50 ?25 
+ . 089^EXP (-199. 999e>;<T ) 739982 v 1 1 B 
TC ( 3 y 1 ) =EXP ( - 1. 99 . 9993 ('T ) 739982 . 1 1 S+EXP ( - . 03322S^ f ) / 4 * 50 723 
-EXP ( - * 05576:1<T >74. 506 V4 

rC (4 y 1 ) =EXP ( "5 .J!cT ) 74788 . 6+EXF ( - . 03322>f{T> 722. 3365 

-EXP ( -199 . 9998^7 ) 77796506 . 55-EXP ( - . ()5576>{;T > /22 . 282 4 
TC < 1 y 2 > =— . 37037:!{ ( EXP ( -199 . 9998X T' ) 739932 . 1 18+EXP ( - ♦ 03322.1“! >74. 5072 
-EXP ( - . 0557 6t7) /A . 50674 ) 




0^ 
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TC ( 2 y 2 ) =2 ♦ 47492>KEXP ( - ♦ 05S76>}cT ) ~1 * 47448:KEXP ( ~ ♦ 03322>i<T > 

1 -♦00044S25i«EXP(“:l.99»999S>}{T) 

TC ( 3 ? 2 ) --= » 05576>i'<EXP ( - ♦ 05576^1 ) /4 ♦ 50674 
- ♦ 03322«EXP < - ♦ 03322^1 ) /4 . 50725 
-199»9998>i<EXP(-199*9998>icT)/39982*118 
TC ( 4 y 2 ) = , 05576>!<EXP ( - ♦ 05576>!<T ^ /22 ♦ 2824+ ' 

+ 199 ♦ 9998)i{EXP ( "1 99 . 9998.1CT ) /7796506 . 55~5 ♦ »EXP ^ “5 . >KT ) /47S8 ♦ 
- * 03322>i{EXP ( - * 03322iKT > 722 . 3865 
TC ( 1 y 3 ) = ♦ 37037^< ( ♦ 03322>KEXP ( - * 03 3 225KT ) /4 ♦ 50725 
+ 199 » 9998>KEXP (-199 ♦99985KT) 739982 *118 
- » 05576)i<EXP ( - * 05576>KT ) 74 ♦ 50674 ) 

TC(2y3)=-17*7966J}:< *03496#EXP(-,05576>KT)/4»50674 
-♦01 242:KEXP ( ~ ♦ 03322:i<T > 74 , 50725 
-199 ♦ 9998JKEXP (-199 ♦9998>f:T) 739982. 118) 

TC ( 3 y 3 ) =1 ♦ 000445>lcEXP ( -199 . 9998>!cT ) + . 00024484>KEXP ( - . 03322^1 ) 

- . 0006899>i{EXP ( - . 05576^+ ) 

TC < 4 y 3 ) =25 . HcEXP ( -5 . )KT ) 74788 . 6+ ♦ 001 lOAtEXP ( - ♦ 03322tT ) 722 . 3865 
- . 0031 liKEXP ( - ♦ 05576JKT ) 722 ♦ 2824 
-39999 . 925KEXP ( -199 . 9998:{cT) 77796506 .55 
TC(ly4)=0.0 
TC(2y4)=0.0 
TC(3y4)=0.0 
TC(4y4)=EXP(-5.)KT) 

RETURN 
END 

STOP — 
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Table A. 7 Otolith state equations 


SENSOR UPDATE 


5: = Ax + Bf 
y ~ Cx + SFR + n 


where 


A 




“ “ 

0 0 

B = 

0 

-40. -200.2 


1 

f 


• -* 


$ 


= [si - a1 ' 


-1 ^ 


(s+200) (s-t-,2) 


C = [1800 ISOOO] 


s+200. 2 
40 


-1 


KALMAN FILTER UPDATE 


k = Ax + Bf 
■y = Cx + SFR + n 

where 



0 10 


o' 

A ■ 

-40. -200.2 1 

B = 

0 


0 0-1 


1_ 

C = 

[l 800 18000 o] 





(s+1) (s+200) (S+.2) 


(s+1) (s+200. 2) -(s+1) 

40 (s+1) s(s+l) 

0 0 


-s 

(s+200) (S+.2) 


S'®®” 
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Table A, 8 Canal state equations 


SENSOR UPDATE 


it = Ax + Bw 
y = Cx + SPR + n 

Where 



0 

1 

0 


0 

A = 

0 

0 

1 

B = 

0 


-.37037 

-.17.7966 

-200.0888^ 


1 

C = 

-23.5785 

-1131.89 

-6371,86] 




$’= [si - aJ"^ = 

(s+199.9998) (s+. 03322) (s+. 05576) 


(s+199.9998) (s+. 08898) 
-.37037 
.37037s 


-(s+200.0888) 1 

s(s+200.0888) -s 

-17.7966(s+.0208) s^ 


KALMAN FILTER UPDATE 


where 


A = 


i = Ax + Boj 


y 

= Cx + SFR 

+ n 
1 




A. 


* 




0 

1 

0 

0 


0 

0 

0 

1 

0 

B =• 

0 

-.37037 

-17.7966 

-200.0888 

1 


0 

0 

0 

0 

-5 

/ ■, J 


1 



Table A. 8 concluded 


C = [-23.5785 -1131.89 -6371.86 63.662o] 

4 = [si - A]“^ = ^ 

(s+5) (s+199.9998) (s+. 03322) (s+. 05576) 


(s+5) (s+199.9998) (s+. 08898) 
-,37037(s+5) 
.37037s(s+5) 

0 


-(s+200, 0888) (s+5) (s+5) 

s(s+5) (s+200,0888) -s(s+5) 

-17.7966(3+5) (s+. 0208) s^(s+5) 


-1 



0 


(s+199.9998) (s+. 03322) 
(s+. 05576) 


qQ 

S| 

§1 
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A«5 Kalman gains subroutine 

The next listing is a program which calculates Kalman 
gains for the otolith system. This program calculates only 
the utricle gains - remember that the saccule gains are twice 
those of the utricle. This iterative routine makes the fol- 
lowing calculations until the Kalman gains reach a steady state: 

1. Calculate state estimate. 

2. Calculate propagated error covariance, 

3. Update state estimate, 

4. Update error covariance, 

5. Calculate error covariance. 

Once steady state gains are obtained (in this case, after 
240 iterations) , they must be tested in the main program for a 
known input-output case,. This is for the purpose of tuning 
the model. The following iterative procedure is used: 

1, Run main program with Ormsby update intervals of .1 
and 1,0 seconds for a known input-output case ( for 
example, a constant yaw acceleration of 1.5°/second 
for 120 seconds, then a sustained yaw rate of 180°/ 

second for 120 seconds) , Plot WTOT for this case. 

' f f .-r. 

2, ' Calculate new Kalman gains for the desired -intervals, 

3, Run main program with same stimulus for new gains. 

Plot WTOT. 


4. If the plots do not match, vary the input variance. 
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measurement noise variance and input power (QU, QM and 
D) to calculate new gains. 

5. Continue this process until the WTOT plots are similar. 

Note that the gains are changed for the otoliths only. This 
was done for simplicity, since the Kalman filters for the canals 
do not change the afferent firing rates appreciably. 

Table A. 9 lists the variables used in the Kalman gains 
routine and their definitions. 


ORIGINAL 
OF POOR QUALEP^i 
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CALCULATES KALMAN GAINS FOR OTOLITH SYSTEM 

DIMENSION PX<3s-3) S-PPC3S-3) ,TM(3s>3) i-TH(9) s-C<3) yGK(3) vS(3/3) i-TT(9) 
DATA C /0»1SOOOOOOE+0470»18000000E+05!>0«OOOOOOOOE+<)0/ 

DATA PX /O*lOOO0OO0E”OlrO*0OO0OO0OE+0070,OO0OOO0OE+00» 

• 1 O»OOOOOOOOE+OO7O*iO0OOOOOE~Ol?O*OO0OOOOOE+OO» 

2 0 ♦ OOOOOOOOE+00 7 0 . OOOOOOOOE+00 7 0 ♦ lOOOOOOOE-01/ 

DATA D /»22500000E+00/ 

DATA QU /0*400/ 

DATA QM /12. 8000000/ 

DATA DT /♦ 250000000/ 

DATA NIT /240/ 

CALL ASSIGN(307 'DKliKALMAN»OTO' ) 

CALCULATE TRANSITION AND COVARIANCE 
MATRICES . 

CALL STM0(DT7TM) 

WRITE(3072) ( (TM(l7 J) 7l=l73) 7j=l?3) 

F0RMAT(' TRANSITION MATRIX ' fZE14*7/19X?ZElA*7/l?XfZElA*7) 

CALL ICM0<DT7TT) 

DO 3 JTH=1?9 

TH ( JTH ) =TT < JTH ) ^QU;:<D)iOK2 

URITE(3074> <TH(I) 7l=l79) 

F0RMAT(^ COVARIANCE MATRIX ' 73E14+7/19X73E14+7/19X73E14»7) 

T=0*0 

DO 85 M=l7NIT 

CALCULATION OF S=PX5{cTM TRANSPOSE 

T=M>KDT 
DO 10 1=173 
DO 10 J=l73 , 

S(l7 J)=0*0 
DO 5 K=ly3 

S(l7J)=S(l7J) +PX ( 1 7 K ) ;}:TM ( J 7 K ) ■ 

CONTINUE 

CALCULATION OF PROPAGATED ERROR COVARIANCE MATRIX 
..i pp=TM5!<S+TH 

IC=0 

DO 20 1=173 
DO 20 J=ly3 
IF(I»GT*J) GO TO 20 
IG=IC+1 

IF<IC*E0*7> GO TO 22 
PP<l7 J)=THaC) 

DO IS K=l73 

PP ( 1 7 J ) =PP < 1 7 J ) -f TM < 1 7 K ) >KS ( K 7 J ) 


15 
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STOP 


20 PP< Ji-D^ppay J) 

c 

C CALCULATION DP Stl ? 1 ) = ( C^PP:i<C+OM ) 

C 

22 S<ii-lJ=^QiT4i-:'.2 

DO 25 I-=l;>3 

S ( I j 2 ) =PP (I y 1. ) >KC ( J. ^ +PP <1 y 2 ^ :i<C < 2 ) +PP < I y 3 ) >KC < 3 ) 
25 Savl>=S<l yi)-fCa):itS(Iy2) 

C 

C CALCULATION OF KALMAN GAINS GK^PPJkC/S ( 1 y 1 ) 

C 

no 30 I=Ly3 

30 GK(I)=Sav2)/SayJ.) 

C 

C CALCULATION OF S<I y J) = (I=GK>^C) 

C 

no 45 1=1 7 3 
no 45 J=ly3 

45 S<IyJ)~"GK<I>»C< J) 

no 50 1 = 1 70 

50 S<l7l)=3<:7 I. >'Fi * 


c 

C CALCULATION QF NEU ERROR CDUAFaANCE MATRIX 

C PX=3.’KPP 

C 

no 4)0 1 = 1 7 0 
no 60 J-i. "3 
■ PX(l7J)=0.0 



no 55 K=; 

55 

PX ( I J > =? 

60 

PX ( J V 1 Jt =r 

C 

Ni? TTE •% 30 

70 

FORMAT ( ' 


1 ■ 


WRIIECSC' 


FORMAT f6F 

85 

CONVINUi-: 

STOP 

ENIi 


4 '/ 


origin Al- 

OP POOR 


pa-GE is 

©jAUsa 
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Table A. 9 Variables used in Kalman gains programs 


C(I) 

D 

DT 

GK(I) 

NIT 

PP(I,J) 
PX(T,J) . 
QM 
QU 

TH(I,J) 

TM(I,J) 

TT(I,J) 


C matrix (see Table A, .7) 

input power 

update interval 

Kalman gains 

number of iterations 

propagated error covariance matrix 

system covariance matrix 

variance of measurement noise 

power of input noise 
2 

TT*QU*D 

state transition matrix 
input covariance matrix 
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A. 6 Kalman gains subrou-tine library 

The last listing is that of the subroutines called by 
the Kalman gain program* Subroutines STMO and STMC were de- 
scribed in section A, 4. Subroutine ICI40 calculates the oto- 
lith input covariance matrix used by the Kalman gains routine, - 
It implements the following equation; 


where 




ii = 


2 T T 
0 B Q B $ 


dt 


= input covariance matrix 

^ = otolith state transition matrix (see Table A. 7) 

B = otolith B matrix (see Table A, 7) 

Q = variance of measurement noise 

The variables used in these subroutines are defined in the 
listing. 


ORIGINAL PAGE IS 
OP POOR QUALITY 
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SUBROUTINE SfilOaiT»TO) 

SUBROUTINE TO CALCULATE OTOLITH STATE 
TRANSITION MATRICES 

DT - UPDATE INTERVAL FOR KALMAN FILTER ' 

TO - 3X3 STM FOR OTOLITH KALMAN FILTER SYSTEM 

DIMENSION T0(3r3) 

CALCULATE TO 
T=DT 

TO ( 1 ? 1) = < 200 » :}CEXP ( “ ♦ 2«T ) - ♦ 2>KEXP ( -200 ♦ 5KT ) ) / 1 99 ♦ 8 
TO ( 2? J. ) = ( EXP ( - , 2>i<T ) -EXP ( -200 ♦ ^;T ) ) /199 ♦ 8 
TO < 3 y 1 ) “ ( EXP ( -200 ♦:^:T ) ) /39760 ♦ 2+ ( EXP ( - ♦ 2'^J ) ) /159 » 84 
1 -(EXP(-T) )/159*2 

T0<ly2>=40>i«<EXP(-200.m-EXP<-*2>KT>)/199»8 
TO(2i-2> = (200wi4EXP(-200»^T)-*2ii(EXP<-*2>KT) )/199»8 
T0(3s-2) = <.£XP(-T) )/159*2-( ♦ 2:KEXP ( - ♦2;{<T ) ) /159 ♦ 84 
1 -<200»>!<EXP(-200«iJ:T) )/39760*2 

T0‘<ly3>--=0a) 

T0<2y3)=0 -0 
T0(373)=EXP(-T> 

RETURN 

END 


SUBROUTINE STAC (DT ? TO 

‘SLiHROUTENE TO CALCULATE SEMICIRCULAR CANAl 

STATE TRANSITION MATRICES 

DT - UPDATE INFEkMaL FOR KALMAN FILTER 
TC - 4X4 STM FOR CANAL KALMAN FILTER SYSTEM 


DIMENSION TCC4p4) 

CALCU! A TE TC 
T--DT 

TC (1 y 1 ) = t . 0C'37 s;LEXP ( - . 03322:;<T ) - , 03322;kEXP < ~ , 05576::< F ) ) / . 02234 
TC ( 2 ? 1 ) =-200 » 033^EXP ( - » 0S57.S:!iT ) /4 *‘S0674 

+ 200 . 005 Se:!<EXP ( - . 03322 :KT > / 4 , 30723 
+ . 039 tEXP ( - 1 99 , 9990 M ) / 399 B 2 » 1 :l 8 
TC,( 3 -y.l )-=EXP ( - 1 . 99 » 9998 ;{!T ) X 39982 . 1 . 18 +EXP ( - ♦ 03322 )icT ) / 4 « 30725 
. -EXP<-> 055765 ;<T )/ 4*30674 

'tC ( 4 y 1 =EXP ( -3 , +'T ) /4 788 . 6+EXP ( ~ . 03322+T ) /22 , 3865 

-.-.XP t - L99 - V998-L f ) / 7796506 ♦ 35-EXF' < - ♦ 03576in' ) /22 » ..'824 
TC ( 1 y 2 > =- .37037;:': f RXP’ < -1 99 * 99V8'+ f > /399B2 ♦ 1 18+EXP < - . 0332240f ) /4 ♦ 3072 
-EXP ( - . 03576;}: D /4 . 50674 ) 

TC ( 2 y 2 ) =2 , 47492;KFXP < - . 0S576;k F) -1 . 47448;icEXP < - ♦ 03322:'i<T ) 

- . L'OO 4 'VolCFEXP ( - 1 99 . 9998;}; f ) 

T C ( 3 y 2 > = . 0 r .5 7 A E .X r-' ( - , O 9. 5 7 6 ;i< T > / 4 . 5 0 6 7 4 
- . 0.'.^227EXP < - , c}3322;LT ) 74 . 30 ’23 


1 
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2 “199 . 999PSEXP ( - J. 99 . 999B;KT ) /399B2 » IIB 

TC (4 » 2 ) = * OSSTd-'KHXP < “ > Q!n576;i<T ) /22 » 2824+ 

1 +199 . 999B;i'l£XP ( - 1 99 . 9998^1 ) /7796506 ♦ , >KEXP ( -5 . 5^:T ) 74788 + 6 

2 “ + 03322JKEXP ( “ ♦ 03322?i<T ) 722 ♦ 3865 

TC < 1 » 3 ) =-■+ 37037;:- C > 03322>i'F:XP ( ~ + 03322>f:T ) 74 » 50725 

1 +199 + 9998;KEXP ( - !.99 + 999B>KT ) 739 9B2 + J. IB 

2 - + 05576;KEXP ( - ♦ 05576>KT ) 74 . 50674 ) 

TC < 2 !- 3 ) =“17 + 7966-K ( - 03496;KEXP < ~ . 05576>KT ) 74 , 50674 

1 - ♦ 0 1 24 25KEXP < - * 03322>!<T ) 74 ♦ 50725 

2 -199 . 999S;iiEXP (-199. 999B:KT ) 739982 + 113) 

TC < 3 j 3 ) =1 + 000 45;KEXP ( -199 + 999Qi^T ) + » 00024484;SEXP C - + 03322>KT ) 

1 “ + 0006899’-{{£XP<“ + 055765KT) 

TC < 4 ? 3 ) =25 + ^l^EXP ( -5 + rT ) 74788 + 6+ + 00 J 104>i{EXP < - + 03322:i<T ) 722 + 3865 

1 - * 0031 l>i<EXP < - * 05576>KT ) 722+2824 

2 -39999 * 92>1:EXP ( -199 + 999S;KT )/7796506 + 55 
TC(l»4)=0+0 

TC(2?4)=0+0 
TC(3p4)=0v0 • 

TC<4r4)=EXP(-5 + ;i<T) 

RETURN 

END 


SUBROUTINE ICHO(DTj-TT) 


SUBROUTINE fU LhLCULATE OTOLITH INPUT 
CQOARImrCE iiuTRtX 

DT - UPDhTE interval FOR KALHAN FILTERS 
Tf “ INPUT CCiiA'rR uANCE HAFRIX 


1 

I 

I 

1 

. 1 . 


DIMENSION rf<-9> 

T=DT 

TT ( 1 ) =“EXP ( - + 4XT 'i ;L9* 7351953E-05 + EXP ( -1 * 2:!^ ) ;{(6 > 5496806E-05 
-EXPC-2 0-rr>;n +y728037E-05 + 5 + 20S3184E-05 
TT(2) = -EXP<-i ,2-KD>:;3 P”:'^-8Q89E-05 + EXP(- + 4:KT .--{(i + 9S70367E-05 
+EXP(-2 + v;KT)^I + 9'/'2B037E-05 - + 000031 SE-OS 
T r ( 3 > =-EXI-- ( - i * ? r r - 2 i3546!;E“03 + EXP ( -2 v 0>KT ) :}-'3 + 1 407035F-03 

+2 + 0728 4 3-'or—>i 


TT ' 4 ) 


•EXPC-i^O 


, ■? 7 2 P 03 7 E - 0 5 + E X P ( - 1 + 2>K T ) +: 1 + 3 0 99 3 6 3 L‘ - 0 5 


" EaP !>.- t 4:; f ) ',.3 , 9 i 40 735E-06 
TT(5)=-E-XPv-2 0;:<r)73 > l40;'0.i5L“03’ + EXP ( -I ♦ 2>KT ) TI + 0 427093!: -03 
+2 + 0<J799'rIL-0A- 
rTC*)=— EX! <-2.0^ i ) + ,3 j ,5 


TT < 7 ) -=0 * 0 
TT(P' -OvO 

TT(9)=0.0 
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STOP 


RETURN 

END 
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